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CHAPTER  1 


INTRODUCTION 


Since  the  early  1970* s,  tremendous  growth  has  been  seen  in  the  development  of 
computer  software  for  weapon  systems.  Part  of  this  development  is  a  result  of  the 
development  of  microprocessors  and  distributed  processing  and  networking.  With 
miniaturization  of  components  in  computer  systems  came  the  ever-increasing  role 
of  software  in  the  weapon  systems .  The  software  is  called  upon  to  perform  more 
complex  tasks  than  ever  before,  e.g.,  weapon's  coordination,  scheduling,  and 
control.  This  increased  role  has  also  meant  a  dramatic  increase  in  software 
costs.  In  1975,  it  was  estimated  that  software  costs  exceeded  hardware  costs  by  a 
factor  of  three  or  four  for  U.S.  Air  Force  Weapons  systems.1  In  1977,  the  costs 
of  software  alone  to  the  entire  U.S.  economy  ranged  from  $10  to  $19  billion.2 

With  such  an  increasing  reliance  on  computer  systems,  there  is  a  major  prob¬ 
lem  in  developing  "error-free"  programs.  For  large  scale  real-time  embedded  com¬ 
puter  systems  such  as  the  TRIDENT-I  Fire  Control  System  (TFCS)  and  its  follow-on, 
TRIDENT-II,  it  is  an  impossible  task  to  check  every  conceivable  logic  path  in  the 
computer  code  for  every  combination  of  possible  inputs,  to  discover  "programming 
errors."  Researchers  and  practitioners  of  software  code  development  have  looked 
to  various  tools  to  cut  down  on  the  number  of  errors  in  the  design  and  develop¬ 
ment  stage.  Included  among  these  tools  or  approaches  are:  structured  code,  a 
"top-down"  approach  to  the  software  design,  and  the  development  of  a  number  of 
automated  verification  and  validation  (V&V)  tools  for  program  checkout;  however, 
they  have  not  proven  to  be  a  complete  answer.  A  quote  from  the  July  1973  Air 
Force  Magazine  states: 

"The  world's  most  carefully  planned  and  generously  funded 
software  program  was  that  developed  for  the  Apollo  series  of 
lunar  flights.  The  effort  attracted  some  of  the  nation's 
best  computer  programmers  and  involved  two  competing  teams. 

Checking  the  software  as  thorough  as  the  experts  knew  how 
to  make  it.  In  the  aggregate,  about  $600  million  was  spent 
on  software  for  the  Apollo  program.  Yet  almost  every  major 
fault  of  the  Apollo  program,  from  false  alarms  to  actual 
mishaps,  was  the  direct  result  of  errors  in  computer  soft¬ 
ware."3 

Another  problem  that  the  U.S.  Government  and  the  Department  of  Defense  (DOD) , 
in  particular,  are  facing  with  software  procurement  is  the  inability  to  establish 
and  enforce  software  reliability  goals  from  contractors.  "How  does  one  develop 
guidelines  or  standards  that  can  be  used  to  ensure  a  certain  quality  in  the  soft¬ 
ware  as  is  currently  imposed  in  military  standards  for  hardware  development?" 
A  more  basic  question  that  needs  to  be  answered  first  is:  "What  is  meant  by  soft¬ 
ware  reliability?"  First,  specific  goals  or  objectives  need  to  be  established. 


This  report  defines  software  reliability  as  "the  probability  that  a  given 
software  program  will  operate  without  failure  for  a  specified  time  in  a  specified 
environment."  The  specified  environment  is  particularly  emphasized  as  it  consti¬ 
tutes  one  of  the  major  assumptions  for  many  of  the  reliability  models  discussed  in 
this  report.  If  the  testing  environment  is  quite  a  bit  different  from  the  actual 
operating  environment,  the  program's  reliability  cannot  be  accounted  for  in  that 
environment.  Software  error  or  failure  is  defined  as  "any  occurrence  attributable 
to  software  in  which  the  system  did  not  meet  its  performance  requirements."  These 
definitions  are  consistent  with  the  majority  of  such  definitions  found  in  the 
literature. 

Knowing  the  current  status  of  the  program  reliability  can  determine  when 
testing  should  be  completed  and  the  program  released  for  operational  use.  It  can 
also  aid  the  software  manager  in  determining  how  best  to  allocate  his  limited 
resources  (manpower,  computer  time)  among  the  various  program  modules  for  testing. 
The  current  program  reliability  can  be  used  in  making  decisions  regarding  design 
tradeoffs  between  reliability,  costs,  performance,  and  schedule.  Another  use  is 
in  evaluating  various  software  engineering  approaches  or  tools  to  find  the  one 
that  leads  to  the  "most  reliable"  program  with  (hopefully)  the  minimum  cost.  The 
literature  is  sadly  lacking  on  controlled  studies  which  indicate  the  performance 
of  software  tools  in  eliminating  errors  in  software  code. 

The  purpose  of  this  report  is  to  provide  a  survey  of  the  various  approaches 
that  have  appeared  in  the  literature  concerning  the  estimation  or  modeling  of  a 
program's  reliability.  This  report  describes  the  underlying  assumptions  for  each 
of  the  models  and  provides  a  data  requirements  list  for  implementation.  The 
various  models  are  contrasted  with  each  other  and  the  relative  merits  or  drawbacks 
are  also  highlighted.  This  report  provides  a  practical  guide  for  the  implements- 
tiou  of  these  procedures  on  a  software  program.  Finally,  the  report  gives  any 
results  of  studies  undertaken  to  analyze  the  performance  of  these  procedures. 
Unfortunately,  this  is  one  of  the  areas  in  which  little  has  been  done.  Most  of 
these  studies  are  either  based  upon  simulated  data  or  data  sets  for  which  the  data 
were  collected  for  purposes  other  than  reliability  modeling.  As  a  result,  some  of 
the  k»y  assumptions  upon  which  these  models  or  approaches  rest  are  violated.  An 
additional  purpose  of  this  report  is  to  provide  the  assumptions  and  data  require¬ 
ments  for  the  various  models.  Steps  will  be  taken  in  the  software  development  for 
TRIDENT-II  to  ensure  the  compatibility  of  the  data  with  the  model  assumptions. 

Over  the  last  15  years,  these  models  and  estimation  procedures  have  evolved. 
There  are  basically  three  different  approaches  that  have  been  identified  in  the 
literature:  Error  Seeding/Tagging  Models,  the  Data  Domain  Approach,  and  the  Time 
Domain  Approach.  Chapter  2  of  this  report  describes  the  Error  Seeding/Tagging 
Models,  Chapter  3  describes  the  Data  Domain  Approach,  and  Chapter  4  describes 
the  Time  Domain  Modeling  efforts.  Chapter  5  describes  any  studies  and  their  con¬ 
clusions  in  comparing  the  performance  of  these  various  approaches  on  actual  data 
sets.  Finally,  Chapter  6  presents  a  number  of  "quick"  estimates  of  reliability. 

Before  beginning  the  description  of  these  various  approaches,  it  must  be  kept 
in  mind  throughout  this  report  that  software  reliability  modeling  is  just  one  of 
many  tools.  It  cannot  provide  all  of  the  answers  that  the  software  managers  must 


HI  w"'h  Tl  'W"*CW  v  1 ’j  •'»  rK'^_. w  ,  ***  v*- v*1, ,’»'  »^v 


NSWC  TR  82-171 


face.  It  must  be  taken  as  a  bit  of  information,  which  along  with  others,  is  help- 
ful  in  making  a  realistic  judgement  concerning  a  program's  status.  Because  of  the 
current  controversy  about  which  of  the  models  is  best  and  because  of  the  uncer¬ 
tainty  about  the  performance  of  the  software  reliability  modeling  approaches,  it 
is  emphasized  that  the  model  that  is  best  suited  to  the  data  be  applied.  The 
resulting  estimate  of  reliability  may  be  used  as  another  source  of  information  in 
determining  program  status . 
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CHAPTER  2 

ERROR  SEEDING/TAGGING  MODELS 


This  approach,  first  proposed  by  Mills, *  involves  "seeding1  a  given  program 
with  a  number  of  known  errors.  The  assumption  is  made  that  the  distribution  of 
the  "seeded"  errors  is  the  same  as  the  distribution  for  the  inherent  errors  in 
the  program.  The  program  is  then  given  over  to  a  testing  team  for  V&V.  Some  of 
the  errors  discovered  by  the  testing  team  are  seeded  errors  while  others  are  in¬ 
herent  in  the  program.  Using  these  counts,  the  total  number  of  errors  inherent  in 
the  program  can  be  estimated.  In  particular,  if  there  are  N  errors  inherent  in 
the  program  and  n  are  randomly  inserted  with  r  errors  being  subsequently  detected 
by  the  quality  assurance  (QA)  team  [k  (k<r)  being  seeded  errors],  it  can  be  shown 
that  the  maximum  likelihood  estimate  (MLEj  of  N  is: 

K  =  [”<r  '  k)1  .  C2-1) 

with  [  ]  being  the  greatest  integer  function. 

The  biggest  drawback  to  this  Seeding  Approach  is  the  assumption  that  is  made 
about  the  distribution  of  seeded  errors  being  the  same  as  the  distribution  of 
inherent  errors.  This  is  an  impossible  assumption  to  check,  especially  in  the 
latter  stages  of  program  development.  At  that  point,  many  of  the  easy  errors 
(e.g.,  misspelled  output)  have  been  eliminated  and  the  only  remaining  errors  are 
the  very  subtle  errors  which  are  extremely  difficult  to  uncover. 

Another  approach,  proposed  by  Rudner,5  avoids  this  problem  by  employing  a 
"two-stage"  or  "two-team"  testing  procedure.  The  program  is  first  given  to  one 
team  for  testing  which  finds  n  errors.  The  program  is  then  turned  over  to  a 
second  testing  team  which  discovers  a  total  of  r  errors,  (k  of  which  were  also 
found  by  the  first  team.)  Using  the  hypergeometric  distribution,  the  MLE  for  the 
total  number  of  errors  in  the  program,  N,  can  be  shown  as: 


where  again  [  ]  denotes  the  greatest  integer  function. 

In  an  article  by  Schick  and  Wolverton,6  reference  is  made  to  a  pair  of  papers 
by  Basin7’8  in  which  the  following  approach  is  taken.  Suppose  a  program  consists 
of  M  statements  from  which  n  are  randomly  selected  and  errors  are  introduced.  It 
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r  statements  are  then  randomly  chosen  and  tested  with  ki  having  inherent  errors 
and  k2  having  seeded  errors,  then  the  MLE  of  N,  the  total  number  of  errors,  is: 

A 

» -  [k.  -  «•» 

with  [  ]  being  the  greatest  integer  function. 

All  of  these  procedures  stem  from  "capture/recapture"  estimation  techniques 
which  estimate  the  total  number  of  animals  of  a  given  species.  A  "tagged"  set  of 
animals  is  released  into  the  environment  and  after  allowing  the  animals  sufficient 
time  to  disperse,  a  second  capture  is  made.  Based  upon  the  number  of  tagged 
animals  released  and  recaptured,  estimates  of  the  total  population  size  can  be 
made. 


In  applying  these  estimation  procedures,  Schick  and  Wolverton6  warn  that, 
based  upon  preliminary  calculations,  the  tag  ratio  (the  average  number  of  tagged 
errors  in  the  sample)  should  be  greater  than  20  to  ensure  the  estimates  are  close 
to  being  unbiased.  These  estimation  procedures  can  be  applied  at  any  point  in  the 
life  cycle  development  of  a  program  to  estimate  the  current  error  content.  The 
biggest  drawbacks  are  in  seeding  the  errors  and  in  the  employment  of  limited 
resources  in  a  two-team  approach.  Few  organizations  can  afford  the  luxury  of 
duplicate  testing  teams  for  a  given  program  or  even  program  modules.  Generally, 
in  the  life  cycle  development  of  a  program,  if  schedules  start  to  slip,  the  time 
is  made  up  at  the  expense  of  the  V&V  effort.  As  a  result,  when  the  program 
reaches  the  testing  team,  all  available  resources  are  spent  to  quickly  perform 
the  testing  tasks  and  release  the  program  to  the  operational  user.  These  pro¬ 
cedures  f.lso  do  not  provide  time-dependent  reliability  measures  of  the  software, 
which  may  or  may  not  be  a  drawback.  This  is  discussed  in  Chapter  4,  the  Time 
Domain  Approach. 
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CHAPTER  3 

DATA  DOMAIN  APPROACH 


The  Data  Domain  Approach  includes  those  procedures  that  estimate  a  program's 
current  reliability  based  strictly  on  the  number  of  successful  runs  observed 
compared  to  the  total  number  of  runs  made.  Included  within  this  category  are 
procedures  that  try  to  employ  test  inputs  for  the  program  that  are  chosen  accord¬ 
ing  to  probability  distributions  of  anticipated  operational  usage.  The  various 
inputs  to  a  program  are  broken  up  into  categories,  and  probabilities  are  then 
assigned  to  those  categories  which  represent  anticipated  uses. 

For  example,  range  might  be  an  input.  It  can  be  broken  up  into  the  cate¬ 
gories  [0,  1500  nautical  miles  (nmi)],  [1501  nmi,  2500  nmi] ,  [2501  nmi,  3500  nmi] 
and  [3501  nmi  or  more].  Probabilities  are  then  assigned  to  each  category,  based 
upon  the  anticipated  operational  usage.  If  it  is  anticipated  that  about  one- 
fourth  of  all  the  ranges  are  3500  nmi  or  more,  that  category  is  assigned  the 
probability  The  inputs  are  randomly  selected  according  to  their  probability 
distributions  and  the  resulting  test  cases  are  run.  The  estimated  reliability 
is  then  simply  the  total  number  of  successful  runs  over  the  total  number  of  test 
runs. 

The  Data  Domain  procedures  try  to  divest  themselves  of  time  between  error 
occurrence  that  the  models  of  the  Time  Domain  Approach  may  employ.  If  time  is  a 
factor  in  some  of  these  models,  it  represents  the  total  elapsed  time  (either  wall 
clock  time  or  CPU  time)  for  a  testing  session  and  not  the  times  of  error  occur¬ 
rence. 

If  a  random  selection  of  inputs,  which  reflect  the  anticipated  operational 
use,  is  made  and  N  runs  are  made,  with  S  being  successful,  then  the  estimate  of 
the  current  program  reliability  is: 


R  =  S/N.  (3.1) 

Using  this  basic  Binomial  Experiment  Approach,  a  number  of  researchers  have  pro 
posed  modifications  of  this  estimator.  Hecht9  proposed  the  estimators: 

A 

Ri  =  (S/ (NxL) )  (3.2) 

and 


R2  =  (S/ (NxLxW)) 


(3.3) 
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where  L  is  the  number  of  machine  instructions  submitted  and  W  is  the  average 
number  of  bits-per-instruction.  The  modified  estimators  allow  for  differences 
in  exposure  to  failure  (by  normalizing  the  estimator  by  the  program  length)  and 
for  differences  in  programs  operating  on  different  machines  with  different  word 
sizes.  Suppose  the  reliability  estimates  between  two  programs  are  compared.  If 
one  program  is  on  a  large  main-frame  computer  that  takes  a  large  amount  of  time 
to  process  and  the  other  program  is  on  a  small  microprocessor  which  executes 
quickly,  the  results  may  be  misleading  if  the  reliability  is  calculated  using 
equation  (3.1). 


Brown  and  Lipow10  have  suggested  a  modification  of  the  basic  estimator  to 
allow  for  the  fact  that  many  times  in  testing,  the  input  is  chosen  to  "stress"  the 
software  program.  The  resulting  estimate  of  reliability  then  tends  to  be  on  the 
"pessimistic  side."  Their  procedure  is  to  take  the  input  space  and  divide  it  into 
"homogeneous"  regions,  P^,  i=l,...,K.  They  are  homogeneous  in  the  sense  of  fault 


generation.  Suppose  N.  runs  are  made  from  the  partition  region  P.,  and  F.  are 

J  J  J 


F  . 

failures.  The  estimate  of  the  "unreliability"  of  that  region  is  If  the 

j 

probability  in  an  operational  environment  of  drawing  points  from  P.  is  P{P.},  the 

J  J 

unreliability  for  the  entire  input  space  can  be  estimated  as  the  sum  over  all 


regions  of  the  corresponding  unreliability  of  that  region  times  the  probability  of 


drawing  an  input  point  from  that  region,  i.e., 


F 

Estimate  of  the  unreliability  for  the  program  -  V  i  PfP.}.  (3.4) 

i=l  N.  1 

l 

The  estimate  of  the  reliability  of  the  program  is  then  given  as  1  minus  the  un¬ 
reliability.  The  main  drawbacks  of  this  approach  are  the  construction  of  parti¬ 
tion  sets  which  are  homogeneous  with  respect  no  error  generation  and  assigning 
a  probability  that  an  input  point  be  drawn  from  a  given  partition  region.  The 
former  is  impossible  to  determine  while  the  latter  introduces  a  lot  of  errors  in 
the  estimate  based  upon  subjective  judgement. 

Corcoran,  Weingarten,  and  Zehna11  proposed  a  model  which  is  more  suited  to 
hardware  reliability  applications,  but  because  of  its  easy  extension  to  software 
reliability  modeling  and  the  fact  that  it  is  a  generalization  of  the  previous 
binomial,  it  is  mentioned  here.  Suppose  there  are  M  sources  or  types  of  software 
errors  that  can  occur.  And  suppose  a.  is  the  probability  that  if  the  ith  type  of 

error  is  observed,  it  is  corrected,  i.e.,  the  conditional  probability 
P  {error  corrected | ith  type  observed)  =  a^ 

where 


3-2 


If  N  runs  are  made  and  P\  errors  of  the  ith  type  are  observed,  then  the  estimate 
of  the  reliability  of  the  program  is: 


*  ■  I *  S 


(3.5) 


where  S  is  the  total  number  of  successful  runs  and 


yi  = 


if  F.  >  0 
x 

if  F.  =  0 


(3.6) 


This  estimator  can  be  shown  to  be  asymptotically  unbiased  and  its  variance  goes 

to  zero  for  large  N.  One  drawback  for  this  model  is  that  the  M  types  of  error 

sources  have  to  be  known  beforehand,  and  the  most  serious  drawback  is  knowing  the 

a.  's. 
x 

These  next  models  discussed  try  to  combine  not  only,  the  results  of  a  given 
set  of  runs,  as  in  the  previous  models,  but  they  also  try  to  take  into  account  the 
input  space. 

The  first  model  is  one  by  Nelson.12  The  basic  assumptions  are: 

Assumptions 

(a)  a  program  may  be  defined  as  a  specification  of  a  computable  function  F 
on  a  set  E 


E  =  (E.:i=l,...,N) 

which  is  the  set  of  all  data  input  values  needed  to  execute  the  program. 


(b)  Execution  of  the  program  for  each  input  E^  produces  output  F(E^). 

(c)  Because  of  imperfections  in  the  program,  the  program  actually  specifies 
a  function  F'  which  differs  from  the  intended  function  F. 

(d)  For  some  of  the  E^,  the  actual  output  F'(E^)  is  within  an  acceptable 
tolerance  of  the  intended  output  F(E^);  i.e., 

F'(E,)  -  F(Ei)|  <  Ai.  (3.7) 

But  for  some  E,,  the  actual  output  F'(E.)  is  not  within  acceptable  limits; 
i.e. ,  **  ^ 

F'(E.)  -  F(E .)  >  A.  (3.8) 

J  J  J 

and  an  error  is  said  to  occur.  N  may  be  very  large,  but  it  is  finite,'  owing  to 
the  fact  that  only  a  finite  number  of  different  values  can  fit  into  the  word  size 


h*r  l.  li-  ~+  ■i—  %  /«  7**  m  -  '1*1*  • — •  -  .PJfc.  r  iW>V>iaia'U<TL«.''  vU  .V  «  IUi, 
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of  the  computer.  Now  suppose  that  E  is  the  set  of  all  inputs  producing  errors  on 
a  given  run;  i.e.,  e 


E  =  {E.; 
e  j 


F'(EJ  -  F  (E . )  >  A,.} 


(3.9) 


Suppose  there  are  ng  elements  in  this  set.  Then  the  probability  of  the  program 
executing  correctly  if  an  input  is  randomly  selected  from  E  is: 

n 


R  =  1  - 


N 


(3.10) 


However,  the  usual  situation  is  that  the  points  E^  are  not  chosen  randomly.  They 

are  chosen  according  to  some  operational  requirement  which  can  be  represented  as 
a  probability  distribution  over  E.  For  this  distribution, 


p^  =  P(E^  is  selected}. 

Hence  the  reliability  of  the  program  can  be  expressed  in  terms  of  these  probabili¬ 
ties  as: 


R  r=  T.  PiC1  “  y^) 

i=l  1  1 


where 


(3.11) 


if  E.*Ee 

if  E.eE 
x  e 


(3.12) 


If  n  runs  are  made  and  the  inputs  are  chosen  according  to  the  probability  distri¬ 
bution  over  E,  then  the  probability  of  all  runs  being  successful  is: 


R  =  Rn 
n 


Pi(1 


V 


(3.13) 


Nelson12  expands  his  model  by  allowing  for  the  fact  that  usually  runs  are  not 
made  independently  of  each  other,  i.e.,  one  of  the  input  variables  may  be  chosen 
in  ascending  order  from  run- to-run.  To  allow  for  this,  Nelson  redefines  the  prob¬ 
ability  distribution  over  the  input  space  as: 


{E^  is  selected  on  the  jth  run  of  the  sequence}. 


Hence,  for  the  jth  run,  the  probability  of  a  failure  is: 


(3.14) 


N 

^  =  5,  pijV 


with  y^  as  previously  defined. 


(3.15) 


3-4 
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The  probability  that  there  are  no  failures  in  n  runs  becomes: 
\  -  (1  -  Pi)  (1  -  P2) • • • (l  -  Pn). 


(3.16) 


As  in  the  previous  binomial  type  models,  the  estimated  reliability  based  on 
Nelson's  first  model  is  simply: 


r=‘-5  ' 


(3.17) 


where  f  is  the  total  number  of  failures  and  the  input  points  are  chosen  according 
to  the  probability  distribution  over  E.  To  construct  this  probability  distribu¬ 
tion,  Nelson  suggests  that  the  ranges  of  the  various  input  variables  be  broken  up 
into  subranges.  Probabilities  are  then  assigned  to  these  subranges  based  upon 
anticipated  operational  usage. 

In  the  TRW  report  by  Thayer,13  the  model  by  Nelson  is  again  modified  by  tak¬ 
ing  the  input  space  E  and  partitioning  it  into  disjoint  regions,  R^  i=l,...,k, 


E  =  \J  R.  and  R.Ar.  =  <j». 
i=l  1  1  J 

The  probability  that  a  point  is  randomly  selected  from  R,  can  be  calculated  as: 

—  J 


»  58  E  Pi  ; 

j  E.eR.  1 
J  1  J 


(3.18) 


i.e.,  all  the  operational  probabilities  of  the  input  points  falling  into  region 

R.  are  summed  over.  If  the  region  R.  ?s  further  divided  into  two  sets  R.'  and  R." 
J  J  J  J 


where 


R'.  =  {E.eR.AE  } 
J  1  1  j  e} 


(3.19) 


R'.'  =  {E.eR.  Ae  }, 
J  i  J  = 


(3.20) 


the  set  R’.  is  derived  consisting  of  all  input  points  in  R.,  which  result  in  suc- 
J  J 

cessful  execution  of  the  program.  RV  consists  of  all  input  points  of  R.  yielding 

J  1 

failures. 


The  probability  of  a  point  falling  in  R^  is  therefore: 


PR!  =  £  (1  '  yi)pi  ' 


E.eR. 
1  J 


(3.21) 


'•  V*  ^  ^  •  *  W  *  *  h  •  •  .  •  4  •  *  fc.'-t*”  -  *  V  t  M  -  4  -  H  - 


„  ,s\V 

'  V  ^ 


uf jiUt/-*.  *?■  -i*.  uAl  uLi^M*  iZjjL  uLii  tit 
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with  =  0  if  EieRj,  and  1  otherwise, 
is  similarly  computed  as: 


PR..  =  S  Y,P, 
Rj  E.sR.  1  1 
J  t-  J 


The  probability  of  a  point  falling  in  RV 

J 


(3.22) 


The  overall  reliability  of  the  program  can  then  be  expressed  in  terms  of  these 
probabilities  as: 


R  =  i  -  £  YtPi 
i=i  1  1 


(3.23) 


-  £  £  vi 

T?  ~ r»  *•  * 


j=l  E.sR. 
J  i  J 


(3.24) 


1  “  ?  PRV 
J=1  J 


(3.25) 


K 

S  PR' 


(3.26) 


The  input  space  has  been  stratified  into  regions  in  exactly  the  same  manner  as 
described  by  Brown  and  Lipow.10  Brown  and  Lipow  created  their  strata  based  upon 
creating  regions  which  were  homogeneous  with  respect  to  error  generation,  while 
Nelson12  suggested  a  partitioning  based  upon  logic  paths. 

Using  results  from  basic  sampling  theory  for  a  stratified  population,  if  n^ 

runs  are  made  in  region  R.,  and  f.  are  failures,  then  the  estimate  of  R  is  again 

J  J  x 

simply: 


*  K  f . 

R  =  1  -  E  iT1  PR  . 

j=l  j  j 


(3.27) 


provided  the  input  points  in  R.  are  chosen  according  to  the  probability  distribu¬ 
te 

tion  over  the  input  space  E.  The  variance  of  this  estimator  can  be  shown  to  be 


3-6 
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«  f  .j 

i 


’  P 

.1  k< 


■3  $ 

*i  M .' 


*  K  PR'.  PR'.' 

Var  {R}  =  £  Jn  J  . 
j=l  j 

This  variance  can  be  minimized  by  taking 


(3.28) 


n.  =  n 
J 


PR'  PR'.' 


PR!  PRV 


(3.29) 


As  in  the  previous  Nelson  models,  the  biggest  drawback  is  the  establishment  of  the 
distribution  over  the  input  space  to  determine  Pg,  and  Pg„ . 

j  j 

A  paper  by  Sugiura,  Yamamoto,  and  Shiiuo14  also  considers  the  binomial  model, 
but  views  the  input  space  as  being  composed  of  two  parts.  The  input  space  is 
broken  up  into  a  user  space,  the  space  where  the  input  is  actually  drawn  from,  and 
its  complement  (Figure  3-1). 

E  ^E’^E 


FIGURE  3-1.  THE  INPUT  SPACE  E  =  Eu  (USER  SPACE KAu. 

The  area  is  never  tested  in  the  V&V  stage.  For  the  entire  input  space  E,  sup¬ 
pose  there  are  a  total  of  Ng  input  points,  of  which  a  proportion  Pg  of  them  re¬ 
sult  in  errors.  The  real  reliability  of  the  program  is  therefore: 


Rg  =  1  -  pE  . 


(3.30) 


However  test  inputs  are  obtained  and  sampled  from  the  user  space  Eu,  where  there 
are  a  total  of  Nu  points,  of  which  a  proportion  Pu  result  in  errors. 


i !.  ■5.-,*t£:v**-UT—  ■-.' 
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Now  suppose  a  set  of  n  runs  are  made  and  f  result  in  failures.  P  can  be 
estimated,  as  seen  before,  as 

P  ~  —  for  n  large.  (3.31) 

u  n  6 

Now  suppose  that  m  bugs  are  eliminated  by  debugging.  Then  the  following  expres¬ 
sion  is  derived 


j-m 


P 


u 


N. 


(3.32) 


where  P^  is  the  failure  ‘probability  before  debugging  and  P^  m  is  the  failure 
probability  after  debugging.  The  points  P^  and  P^  m  are  estimated  using  equation 
(3.31)  on  two  separate  teisting  occasions.  Equation  (3.32)  can  be  plotted  as  a 
function  of  m  by  holding  fixed.  Calculating  the  slope  of  this  line,  Nu  can  be 
estimated  as  > 


A 


1 

A - 

Slope 


(3.33) 


Least  squares  can  also  be  used  to  fit  equation  (3.32)  to  data.  The  software  be¬ 
comes  perfect  when 


pj-m  _  pj  _ 
u  “  u 


u 


(3.34) 


i . e . ,  when 

m  =  PJN  ,  (3.35) 

u  u  ’ 

which  simply  is  when  the  number  of  errors  removed  is  the  same  as  the  total  number 
of  input  points  leading  to  failures.  This  assumes  no  new  errors  are  introduced  in 
the  debugging  process.  If  there  are  bugs  which  have  complex  causes,  additional 
software  errors  might  be  introduced  when  correcting  those  errors.  Suppose  that 
such  bugs  exist  in  equal  probabilities  in  the  input  space.  Then  the  probability 
of  such  bugs  residing  in  the  user  space  is  N  /N„.  If  such  bugs  in  their  correc- 

tion  produce  B  new  errors,  the  failure  probability,  after  the  elimination  of  m 
errors,  is: 

pr=tu -rt'-icM-  <3-36) 

u  E 
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The  program  is  error  free  when 


(3.37) 


(3.38) 


where 

N  B 

D*ir-  (3.39) 

If  the  program  is  at  a  particular  point  in  the  testing  stage  where  the  failure 
probability  is  P3 

then 

PJN 

m  =  (3.40) 

is  called  the  "remaining  bug  index."  The  unknowns  in  equation  (3.36)  can  be  esti¬ 
mated,  using  the  results  of  several  testing  sessions  and  equation  (3.31). 

This  is  a  very  simplistic  model  that  has  not  been  employed  on  any  data  sets. 

Their  formulation  also  rests  heavily  on  obtaining  good  estimates  of  the  P3 ' s 

using  equation  (3.31)  in  order  to  fit  the  straight  line  as  a  function  of  m.  This 
means  that  the  number  of  runs  for  a  given  testing  session  should  be  quite  large. 

The  next  reference  in  the  Data  Domain  section  is  to  a  paper  by  Elliot,  et. 
al.15  The  techniques  in  their  paper  again  employ  the  simplistic  assumption  of  a 
binomial  experiment.  For  this  reason,  they  assume  that  the  runs  are  independent 
and  the  true  reliability  is  the  percentage  of  points  in  the  input  space  which 
result  in  the  program  running  correctly,  i.e., 

R  =  |  ,  (3.41) 


3-9 
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where  S  is  the  total  number  of  points  in  the  input  space  E,  which  result  in  suc¬ 
cessful  execution  of  the  program.  N  is  the  total  number  of  points  in  that  space. 

They  propose  two  testing  procedures  to  determine  whether  a  program  has  reach¬ 
ed  a  given  reliability.  One  is  based  on  a  fixed  sample  size,  the  other  is  a  se¬ 
quential  testing  procedure.  The  fixed  sample  size  is  the  usual  hypothesis  testing 

procedure  for  a  binomial  probability.  The  user  specifies  a  size  for  the  Type  I 
error  (the  probability  of  rejecting  the  program  when  it  has  a  desired  reliability 
level)  and  a  Type  II  error  (the  probability  of  accepting  a  program  when  the  reli¬ 
ability  is  no  more  than  a  specified  level).  Using  these  values,  tables  are  given 

which  provide  the  number  of  tests,  n,  to  run  and  the  maximum  number  of  failures, 
f,  that  are  allowed  if  the  program  is  to  be  accepted.  The  program  tester  randomly 
selects  a  subset  of  n  input  points  from  the  input  space  E  and  runs  the  program. 
If  more  than  f  program  failures  are  observed,  the  hypothesis  that  the  program  has 
reached  the  desired  reliability  level  is  rejected;  otherwise,  it  is  accepted. 

For  the  sequential  procedure,  the  tester  specifies:  a  minimum  acceptable 
reliability  Rmin»  a  probability  a  that  the  program  with  this  reliability  will 

pass  testing,  a  probability  Rmax  for  which  one  wants  to  be  "almost  sure"  that  the 

software  will  pass,  and  a  probability  p  that  the  software  with  this  probability 
will  fail  the  test.  The  sequential  procedure  is  to: 


(a)  Accept  the  software  if 

F  <  -h2  +  BNt  ; 

(b)  Reject  the  software  if 


F  >  h!  +  BNt  ; 

(c)  Otherwise,  continue  testing 

where  F  is  the  total  number  of  failures  experienced  up  through  NT  tests  (NT=1, 

2 , . . .  )  , 


hi  = 

[£n(l-a)  -  £np]/D 

> 

(3.42) 

h2  = 

(£n(l-p)  -  £na]/D 

) 

(3.43) 

B  = 

[£nR  -  £nR  .  ]/D 

1  max  min 

y 

(3.44) 

and 

D  = 

£nR  -  £nR  .  -  £n(l 

max  min 

-  R  )  +  £n(l  -  R  .  ) 
max-  ram 

(3.45) 

As  an  example, 
something  less 

by  setting  a  = 
or  larger,  the 

suppose  a  program  is  tested  in  which  if  the 

than  R  .  =  .7  it  is  desirable  to  limit  the 

min 

.05.  On  the  other  hand,  if  the  program  has 
chance  of  rejecting  it  should  be  p  =  .1. 

true  reliability  is 
risk  of  releasing  it 

a  reliability  of  .95 

3-10 


(3.46) 

(3.47) 

(3.48) 

(3.49) 


D  =  2.097 
=  1.074 
h2  =  1.378 

and 

B  =  .146 


(3.50) 

(3.51) 

(3.52) 

(3.53) 


So  the  sequential  procedure  is: 

(a)  Accept  the  software  if: 

F  <  -1.378  +  . 146Nt; 

(b)  Reject  the  software  if: 

F  >  1.074  +  .146Nt; 

(c)  Otherwise,  continue  testing. 

The  advantage  of  the  sequential  procedure  over  the  fixed  sajnple  size  is  that  on 
the  average,  a  sequential  procedure  requires  less  testing  than  a  fixed  sample  size 
to  achieve  the  same  levels  for  the  Type  I  and  Type  II  errors. 

The  last  model  considered  in  this  section  is  LaPadula's  Reliability  Growth 
Model16  (see  also  References  17  and  18).  The  approach  is  to  fit,  using  least 
squares,  a  reliability  curve  through  the  success/failure  counts  observed  at 
various  stages  of  the  software  testing.  More  specifically,  the  assumptions  are: 


Model  Assumptions 

(a)  Testing  is  conducted  in  a  series  of  N  stages.  A  stage  is  marked  by  any 
change  or  modification  to  the  program. 


(b)  At  each  stage,  n^  i=l,...,N  tests  are  performed  of  which  are  success¬ 
ful.  The  number  of  tests  performed  at  a  given  stage  is  not  fixed  in  advance. 
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(c)  After  the  completion  of  the  N-th  stage  (which  itself  is  not  set  in  ad¬ 
vance),  a  growth  curve  of  the  form 


R(k)  =  R(u)  -  A/k 


(3.54) 


is  fitted  to  the  data.  R(k)  is  the  reliability  of  the  program  during  the  kth 
stage  of  testing.  R(u)  is  the  value  of  the  R(k)  as  k-*»  and  A  is  a  growth  pa¬ 
rameter.  If  A  >  0,  the  reliability  of  the  program  increases  while  A  <  0,  the 
reliability  decreases. 

To  estimate  the  two  unknowns  R(u)  and  A,  least  squares  estimates  can  be 
used.  The  desire  is  to  minimize 


S  =  £  R(k)  -  ~ 

k=l\  nk 

N  / 

=  £  (  R(u)  -  i 

k=l  \ 


R(u)  -  A/k  - 

nk 


The  estimates  which  minimize  this  expression  are  found  to  be: 


(3.55) 


(3.56) 


Estimates  -  Least  Squares 


A  =  N 


k=l  nkk 


n  ; 

E  E 

k=l  j 


k  %)(#, ;) 


(3.57) 


n(£  \ 

\k=l  k2 


T  r  -  N  1  N  S, -1 

•  i[‘  5,  i  •  =;] 


(3.58) 


The  only  data  then  required  to  estimate  the  reliability  curve  are: 

Data  Requirement 

The  number  of  tests,  n.,  performed  at  each  stage  and  the  number  of  successes 
observed  at  that  stags,  s^. 

The  relationship  of  the  reliability  of  the  software  to  the  stage  number  of 
the  testing  sequence  is  hard  to  justify.  Moreover,  a  stage  can  have  an  arbitrary 
number  of  tests  composing  it.  The  only  thing  that  marks  the  end  of  one  stage  and 
the  beginning  of  another  is  some  change  to  the  program. 
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CHAPTER  4 

TIME  DOMAIN  APPROACH 


The  Time  Domain  Approach  to  software  reliability  modeling  has  received  the 
greatest  emphasis  in  the  applicable  literature  as  it  does  in  this  report.  This 
approach  attempts  to  utilize  either  the  times  of  error  occurrences  and  the  result¬ 
ing  times  between  error  occurrences  or  the  number  of  error  occurrences  per  time 
period  to  model  the  error  generation  processes.  In  general,  the  models  can  be 
used  to  predict  the  expected  time  until  the  next  error  occurrence  or  the  expected 
number  of  errors  in  the  next  interval  of  testing.  These  models  were  originally 
motivated  by  hardware  reliability  concepts  and  many  of  the  terms  used  in  hardware 
reliability  modeling  are  carried  over  into  the  software. 

Over  the  last  10  years,  many  models  were  proposed  and. extensions  to  them  were 
given.  There  is  still  quite  a  lot  of  controversy  about  which  is  the  "best"  model 
to  use  on  a  software  data  set.  Some  studies  were  done  comparing  the  various 
models  on  simulated  and  real  data  sets  (see  Chapter  5)  and  some  studies  are 
currently  under  way,  but  more  research  is  needed.  The  best  advise  for  applying 
these  models  to  a  software  error  data  set  is  to  apply  a  number  of  them  to  see 
which  appears  to  best  model  the  data;  that  is  one  purpose  of  this  report.  By 
providing  a  general  overview  of  the  various  models,  their  assumptions,  and  data 
requirements;  a  number  of  models,  which  seem  to  be  close  to  the  actual  way  the 
data  was  generated,  can  be  chosen.  By  applying  some  of  these  candidate  models, 
the  best  model  for  a  set  of  data  can  be  established. 

Section  4.1  discusses  some  of  the  hardware  reliability  concepts  and  terms 
that  were  adapted  to  software  modeling.  The  difference  in  hardware  versus  soft¬ 
ware  modeling  is  pointed  out.  Section  4.2  begins  the  discussion  of  software 
modeling  with  some  of  the  classical  adaptations  of  hardware  concepts  to  software 
models.  They  are  classical  in  the  sense  that  many  of  these  models  are  based  on 
an  exponential  distribution  for  the  time  between  error  occurrence  and  the  rate  of 
error  occurrence.  The  latter  is  determined  by  the  number  of  errors  in  the  pro¬ 
gram  at  the  time  of  the  test.  Section  4.3  discusses  the  "Bayesian"  philosophy 
applied  to  software  modeling.  This  is  followed  by  Section  4.4,  which  deals  with 
attempting  to  model  the  behavior  of  the  program  as  a  Markov  process. 


The  models  chosen  in  this  report  were  selected  to  provide  the  reader  with  an 
idea  of  the  numerous  approaches  that  have  been  proposed  for  software  modeling. 
An  extensive  reference/bibliography  is  provided  at  the  end  of  this  report  which 
may  be  of  benefit  to  researchers  in  this  area.  An  excellent  report,  giving  an 
overview  of  software  modeling  in  general  and  containing  an  extensive  reference 
list,  is  a  report  written  by  Gephart,  et.  al.18  This  report  is  highly  recommended 
for  a  researcher  in  this  field. 
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4.1  HARDWARE  VERSUS  SOFTWARE  RELIABILITY  MODELING 

In  hardware  reliability  modeling,19  a  number  of  key  concepts  have  been 
adapted  for  software  modeling.  The  hazard  rate  Z(t)  for  a  component  (software 
program)  is  defined  us  the  conditional  probability  that  a  failure  (error)  happens 
in  an  interval  (t,t+At)  given  that  the  component  (program)  has  not  failed  up  to 
time  t.  If  T  is  the  time  when  a  failure  (error)  occurs,  then: 


Z(t)At  =  P{t<T<t  +  At  |  T>t}  . 


(4.1) 


The  unconditional  probability  provides  the  failure  (error)  probability  density 
function,  f(t),  for  the  component  (software  program);  i.e., 


f (t)At  =  P{t  <  T  <  t  +  At] 


(4.2) 


The  hazard  function  can  be  related  to  the  pdf  of  the  time  of  failure  (error) , 
f(t),  as: 

Z(t)  -  — —  (L  ^ 


(4.3) 


where  F(t)  is  the  cumulative  distribution  of  the  time  to  failure;  i.e., 


F(t)  =  f  f (x)dx. 


(4.4) 


The  function 

R(t)  =  1  -  F(t) 

is  called  the  reliability  function  of  the  component  (program) . 
From  equation  (4.4),  it  can  be  seen  that: 

z(*)d*  =  idr(£-k-v 


(4.5) 


(4.6) 


t  .  t 

f  Z(x)dx  =  -log  [1  -  F(x) ] 


-  /  Z(x)dx  =  log 


(4.7: 


1  -  F(t)  =  R(t)  =  exp  { 


-/ 


Z(x)dx} 


(4.8; 


4-2 


NSWC  TR  82-171 


i 


ft 
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Thus  once  a  hazard  rate  function  for  a  component  (program)  is  specified,  the  reli¬ 
ability.  function  R(t)  is  then  determined.  Once  the  reliability  function  has  been 
established,  the  expected  time  between  failures  (errors)  or  Mean  Time  Before  Fail¬ 
ure  (MTBF)  is  calculated  as: 

oo  oo 

MTBF  =  /  R(x)dx  =  f  tf(t)dt.  (4.9) 
o  o 

Many  of  the  models  in  the  next  paragraph  present  forms  for  the  hazard  rate  which, 
using  the  previous  relationships,  determine  the  reliability  function  for  the  pro¬ 
gram  and  the  MTBF.  Some  of  the  models  in  Section  4.2  contain  terms  that  do  not 
have  counterparts  in  hardware,  e.g.,  the  number  of  errors  remaining  and  the  time 
required  to  discover  the  remaining  errors . 

The  concepts  of  hardware  reliability  modeling  were  adapted  to  software  model¬ 
ing.  This  is  not  to  imply  that  the  behavior  of  software  is  similar  to  hardware; 
quite  the  opposite  is  true.  Software  does  not  wear  out  over  its  life  cycle  as 
hardware  does.  In  the  reproduction  of  software,  there  is  no  generation  of  new 
random  software  errors  introduced  in  subsequent  copies.  Duplicate  software 
programs  yield  identical  results.  Moreover,  software  does  not  change  during 
repeated  operational  use  as  hardware  does.  It  is,  in  fact,  that  inconstant 
property  of  hardware  upon  which  the  probabilistic  modeling  of  hardware  failure 
occurrence  is  based.  For  software,  it  is  the  unchangeability  over  time  that  makes 
software  error  generation  independent  of  time.  The  elapsing  of  a  time  variable 
does  not  cause  software  errors.  For  this  reason,  a  number  of  researchers  have 
strongly  questioned  the  modeling  of  error  occurrence  in  which  time  plays  a  fac¬ 
tor.  (See  References  20  and  21.)  It  is  not  a  direct  relationship  between  a  time 
variable  and  error  generation  that  is  modeled,  however,  but  an  indirect  relation¬ 
ship  as  a  result  of  the  randomization  of  the  input  space  for  a  program  in  opera¬ 
tional  use. 

Within  a  program  are  latent  errors  which  are  discovered  when  a  certain  com¬ 
bination  of  input  variables  cause  execution  of  the  program  to  go  down  the  path  in 
which  the  error  lies.  Because  of  the  very  large  number  of  combinations  of  input 
variables  that  are  possible,  the  operational  usage  of  a  program  gives  the  appear¬ 
ance  of  randomization  over  the  input  space.  This,  in  turn,  causes  the  error 
occurrences  to  take  on  the  appearance  of  following  a  probabilistic  model  over 
time.  It  is  characterization  of  this  probabilistic  nature  on  which  the  modeling 
is  based. 


4.2  CLASSICAL  SOFTWARE  MODELS 


'K’ 


a 


4.2. 1  Weibull  Model 

Since  a  number  of  the  concepts  of  hardware  reliability  theory  were  initially 
adapted  to  software,  one  of  the  earliest  models  to  be  applied  was  the  Weibull 
Model.  Because  of  the  nature  of  the  Weibull  distribution,  it  can  bA  used  to 
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model  increasing,  decreasing,  or  constant  failure  rates  for  software.  The  form  of 
the  hazard  rate  is  taken  as: 


z<«  ■  i  (r 1 


where  a,b  are  constants  >  0  and  t  >  0  (A. 10) 


so  that  if  a  >  1,  the  error  rate  increases  with  time;  if  a  <  1,  it  decreases  with 
time;  and  if  a  =  1,  there  is  a  constant  failure  rate  over  time.  The  corresponding 
pdf  for  the  time  to  failure  is  the  Weibull  distribution;  i.e., 


f(t)  *1(1 


a-1 


exp 


t  >  0 


(A. 11) 


*  J 

|J  kr?v! 


l "  i  *>  :1 

•  *«•»  ■71 

‘*i"  *11 

C.  *\/ 


with  the  cumulative  distribution  function 


I 


£••• 


i 


‘iN 

ft 


F(t)  =  f  f (x)dx  =  1  - 


exp 


(A. 12) 


(Notice  that  if  a  =  1,  i.e.,  a  constant  failure  rate,  f(t)  becomes  the  exponential 
distribution. ) 


The  reliability  function  is  therefore: 

a) 

R(t)  =  1  -  F(t)  =  exp 
and  hence,  the  MTBF  is: 


MTBF 


CO  00 

=  /  R(t)dt  =  /  tf(t)dt  =  I  r(\) 

n  a  '  ' 


(A. 13) 


(A.1A) 


where  F  (•)  is  the  gamma  function. 

Coutinho22  (also  see  Reference  18)  proposes  estimating  the  unknowns  using  as 
input  the  following  data  requirements. 

Data  Requirements 

(a)  The  total  number  of  errors  in  each  time  interval  of  testing 

n£»  i  =  1 , . . . ,K  , 

(b)  The  length  of  the  testing  interval 

d.,  i  =  1,...,K  , 

(c)  The  total  number  of  time  intervals,  K,  and 
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(d)  The  cumulative  number  of  errors  found  to  date 
K 

M=E  n.  . 
i=l  1 

The  estimators  are  obtained  using  graphical  procedures,  method  of  moments,  least 
squares,  or  MLEs.  Tor  instance,  in  the  case  of  least  squares,  let 


m  -  a, 


bQ  =  -a£.n(b) , 


F(i)  =  a 


L*j 

=1  J 


M 


(4.15) 

(4.16) 


(normalized  cumulative  number  of  errors  found  up 
through  the  ith  time  interval),  (4.17) 


and 


x>  =  ■ 

From  the  expression  of  the  cumulative  distribution  function 


we  have 


F(t)  =  1  -  exp  < 


r^tj"  exp 


so 


In 


[r^Fjtj]  =(e) 


and  thus, 


that  is, 


An 


£n 


=  a£n(t)  -  a£n(b); 


y  -  mx  +  b 


(4.18) 


(4.19) 

(4.20) 

(4.21) 

(4.22) 

(4.23) 
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The  estimates  of  the  reliability  function  and  MTBF  are  qiven  as: 


(4.31) 

(4.32) 


Wagoner23  (also  see  Reference  18)  also  applies  the  previous  procedure  to  a  set  of 
software  data,  but  suggests  that  the  d^'s  should  be  measured  in  CPU  time  rather 

than  wall  clock  time.  This  a  good  suggestion  and  should  be  considered  for  all  of 
the  models.  CPU  time  reflects  the  variation  in  testing  effort  from  period  to 
period.  It  also  takes  into  account  when  no  testing  is  going  on.  This  is  dis¬ 
cussed  again  in  relationship  to  Musa's  Model. 


4.2.2  Shopman  Model 

One  of  the  earliest  proposed  software  models  was  derived  by  Martin  Shooman 
(References  24  through  29).  The  basic  assumptions  are: 


Model  Assumptions 


(a)  The  number  of  errors  in  the  code  is  a  fixed  number. 


(b)  No  new  errors  are  introduced  into  the  code  through  the  correction  pro¬ 
cess. 


(c)  The  number  of  machine  instructions  is  essentially  constant  (i.e.,  the 
program  is  relatively  mature). 

(d)  The  detections  of  errors  are  independent. 

(e)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(f)  The  error  detection  rate  is  proport.'  onal  to  the  number  of  errors  remain¬ 
ing  in  the  code. 

Suppose  X  is  the  quantity  of  debugging  time  (in  months)  spent  on  the  system  since 
the  start  of  the  testing  phase  and  suppose  t  is  the  operating  time  (measured  in 
CPU)  of  the  system.  Using  assumption  (f)  at  any  time  t,  the  hazard  rate  is 


where  K  is  the  proportionality  constant  and  £r(x)  is  the  error  rate.  This  is 

taken  as  the  number  of  errors  remaining  in  the  program,  after  X  months  of  debug¬ 
ging,  normalized  with  respect  to  the  total  number  of  instructions  in  the  code. 
This  error  rate,  £r(x),  is  mathematically  expressed  as: 


er(1I)  =  j;  -  «c(t) 


(4.34) 


where  E^,  is  the  total  number  of  errors  initially  in  the  program;  1^  is  the  number 
of  machine  instructions;  and  £c(i)  is  the  cumulative  number  of  errors  fixed  in  the 
interval  from  0  to  x,  normalized  by  the  number  of  machine  instructions.  Since  E^, 
and  1^  are  constant  [assumptions  (a)  and  (c)]  and  since  no  new  errors  are  intro¬ 
duced  in  the  correction  process  [assumption  (b))  as: 


X  «  £  (x)  *  j- 


(4.35) 


er(x)  +  0. 

Combining  equations  (4.33)  and  (4.34), 


Z(t)  =  K 


4  ■ 


Thus,  the  reliability  function  is: 


R(t)  =  expJ-  K 


ft 


and  the  MTBF  is: 


(4.36) 


(4.37) 


(4.38) 


MTBF  = 


(4.39) 


-  ec(i) 


The  only  unknowns  in  this  model  are  E,p  and  K.  These  quantities  can  be  estimated 
in  one  of  two  ways . 

The  simplistic  procedure  is  to  use  the  moment  technique.  The  required  data 
inputs  for  this  estimation  procedure  are  given  in  the  following. 
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Data  Requirements  -  Moment  Technique 

Run  a  functional  test  of  the  program  after  two  different  debugging  times, 
Ti  <  X2,  which  are  chosen  so  that  e  (x i )  <  e  (.X2) ,  and  record  the  following  infor¬ 
mation:  c  c 

(a)  For  each  testing  period,  record  the  number  of  test  runs  that  were  made, 
i.e.,  ri  and  r2  (usually  ri  =  r2). 

(b)  For  each  testing  period  and  for  each  run,  record  the  amount  of  CPU  time 

that  the  program  successfully  executed.  If  out  of  the  r^  runs  made,  m^  were  suc¬ 
cessful  with  execution  times  T\  ,  and  ^  were  unsuccessful,  but  had 

successful  execution  times  of  t.,,...,t.  before  the  errors  were  discovered, 

il’  ’  i,r.-  m.  ’ 

*  1  1 


k- 

:-'V- 

kV\ 


r  .-m. 


H.  =  T.  .  +  E  t.  . 

1  fil  1J  j=l  1J 


(4.40) 


is  the  total  amount  of  successful  execution  time  in  the  ith  functional  testing 
period. 

The  constant  failure  rate  for  the  ith  functional  testing  period  is  then  estimated 
as : 


=  number  of  failures  per  hour 


(4.41) 


r .  -  m. 

i  i 


(4.42) 


Since  the  MTBF  for  a  constant  failure  rate  is  the  reciprocal  of  the  failure  rate, 
the  MTBF  for  the  ith  functional  testing  period  can  be  estimated  as: 


a  H. 

MTBF.  =  - - - 

x  r.  -  m. 


(4.43) 


If  this  expression  is  equated  with  the  expression  for  MTBF,  based  on  the  model,  it 
can  be  seen  that: 


nl  A  i 

— - -  =  MTBF  =  - - - . 

ri"  ml  1  F 

K  -  ec  (t.) 


(4.44) 
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n  r  ; 
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~2-m-  =  mtbf2  = - - - 

r2  m2  rE 

K[^-£o  <*•> 


(4.45) 


m  m 


There  are  two  equations  in  two  unknowns,  so  solving  for  K  and  ET,  the  estimates 
are  obtained. 


Estimates  (Moment  Estimators) 

-  t  [(Z2/Zi)  ec(ti)  -  ec(x2^ 
ET  “  XT  (Z2/Zi)  -  1 


(4.46) 


!*1|  ,*W 


(4.47) 


4 


where 


z  =  !l1 

Zi  H, 


(4.48) 


The  problem  with  this  estimation  procedure  is  the  variation  in  the  estimates 
as  a  function  of  the  two  debugging  times,  Xi  and  X2,  chosen.  Gephart  et.  al.18 
found  that  the  estimates  of  E^,  and  K  varied  quite  a  bit  depending  on  the  two 

chosen  points.  They  suggested  that  a  number  of  pairs  be  chosen  and  the  averages 
of  the  resulting  estimates,  using  the  previous  equations,  be  used.  Using  the 
median  of  these  derived  estimators  as  a  possible  estimate  could  also  be  con¬ 
sidered. 

A  second  estimation  procedure  is  based  on  the  maximum  likelihood  proce¬ 
dure.28*30  The  data  requirements  are  the  same  as  required  for  the  moments  estima¬ 
tion  procedure. 

Data  Requirements  -  Maximum  Likelihood  Estimation 
(Same  as  moments  techniques.) 
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The  MLEs  are: 


Estimates  (Maximum  Likelihood  Estimators) 


[(Z2/Z1)  e  (ti)  -  e  (t2)l 
"  (Z2/Zi)  -  1 


(4.49) 


A 


(4.50) 


Notice  that  E^,  is  the  same  as  the  moments  estimators.  Here  again,  it  is  suggested 

that  a  number  of  pairs  of  points  be  chosen  and  the  average  or  median  of  the  re¬ 
sulting  estimates  of  and  K  be  used. 


4.2.3  Jelinski  and  Moranda  "De-Eutrophication11  Model 

Another  early  model  was  one  proposed  by  Jelinski  and  Moranda31  while  working 
for  the  McDonnell  Douglas  Astronautics  Company.  They  developed  this  model  for 
use  on  the  Navy  NTDS  software  and  for  a  number  of  modules  of  the  Apollo  program. 
As  can  be  seen  in  this  paragraph,  their  work  spawned  quite  a  few  variations  of 
their  basic  model. 

Model  Assumptions 

(a)  The  rate  of  error  detection  is  proportional  to  the  current  error  content 
of  a  program. 

(b)  All  errors  are  equally  likely  to  occur  and  are  independent  of  each 
other. 

(c)  Each  error  is  of  the  same  order  of  severity  as  any  other  error. 

(d)  The  error  rate  remains  constant  over  the  interval  between  error  occur¬ 
rences. 

(e)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(f)  The  errors  are  corrected  instantaneously  without  introduction  of  new 
errors  into  the  program. 


These  assumptions  are  basically  the  same  ones  stated  for  Shooman's  Model. 
(In  fact,  this  report  shows  that  the  two  models  are  equivalent  when  the  correct 
correspondences  are  made.)  The  biggest  questions  are  with  regards  to  assumptions 
(c)  and  (f).  It  is  difficult  to  envision  a  situation  in  which  a  perfect  error 
correction  process  is  achieved.  The  instantaneously  corrected  error  part  of  the 
assumption  can  be  avoided  by  not  counting  errors  which  were  previously  detected, 
but  were  not  corrected.  Assumption  (c)  can  be  avoided  by  dividing  the  errors 
into  classes  based  upon  severity.  For  instance,  one  might  have  a  category  for 
critical  errors,  a  category  for  less  serious  errors,  and  one  for  minor  errors 
(e.g.,  a  misspelled  word  on  an  output).  Software  reliability  models  are  then 
developed  for  each  type.  This  approach  is  suggested. 

Using  assumptions  (a),  (b),  (d) ,  and  (f),  the  hazard  rate  is  defined  as: 

Z(t)  =  <(>  [N  -  (i  -  1)]  (4.51) 

where  t  is  any  time  point  between  the  discovery  of  the  (i  -  l)th  error  and  the  ith 
error.  The  quantity  <p  is  the  proportionality  constant  given  in  assumption  (a). 
N  is  the  total  number  of  errors  initially  in  the  system.  Hence,  if  i  -  1  errors 
have  been  discovered  by  time  t,  there  are  N  -  (i  -  1)  remaining  errors  so  the 
hazard  rate  is  proportional  to  this  remaining  number.  Figure  4-1  is  a  plot  of  the 
hazard  rate  versus  time.  As  can  be  seen,  the  rate  is  reduced  by  the  same  amount 
<J>  at  the  time  of  each  error  detection. 

If  =  t^  -  ti_1,  i.e.,  the  time  between  the  discovery  of  the  ith  and  the 
(i  -  l)st  error  for  i  =  l,...,n  where  to  =  0,  using  assumption  (d),  the  X^'s  are 
assumed  to  have  an  exponential  distribution  with  rate  Z(t^).  That  is: 

f(X.)  =  <))  [N  -  (i  -  1)]  exp  {-<j>[N  -  (i  -  1)]X.}  (4.52) 
z(t) 


N0 


FIGURE  4-1.  DE-EUTROPHICATION  PROCESS 
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so  the  joint  density  for  all  the  X^'s,  using  assumption  (b),  is: 

L(Xlt...,X  )  =  n  f (X. )  =  n  0  [N  -  (i  -  1)]  exp  { -<t> [N  -  (i  -  1)]X.}  .  (4.! 

i=l  i=l 


Taking  the  partial  derivatives  of  AnL  with  respect  to  N  and  <()  and  setting  the  re¬ 
sulting  equations  equal  to  zero,  the  solutions  for  the  following  set  of  equations 
are  obtained  as  MLEs  for  N  and  <j>. 


Estimates  -  Maximum  Likelihood 


4>mt  = 


a  /  n 


*  £  xi  ■  £ (1  -  1)xi 


(4.54) 


i=i  Ha  -  (i  -  1) 


NML“  n1  (]£  (i  -  1)X.\ 

12  X.  Vi=l  7 

i=l  1 


(4.55) 


Equation  (4.55)  is  solved  for  N  using  munerical  techniques  (e.g.,  Newton-Raphson) 
and  is  then  substituted  into  equation  (4-54)  to  obtain  an  estimate  of  <j».  The 
estimate  of  the  MTBF  is  therefore  derived  after  the  jth  error  occurrence  as: 


MTBF  (for  the  (j+l)st  error}  =  ^ —  = 


Z<V  Wm,  ' 


(4.56) 


A  report  by  Tal32  derives  the  least  squares  estimators  for  N  and  <(>  as  the 
estimators  which  minimize  the  sum  of  the  squared  differences  between  the  observed 
time  between  failures  and  their  mean  values,  i.e.,  the  MTBFs.  The  quantity  to  be 
minimized  is: 


t  (X.-MTBF*)** 


(4.57) 


Again  taking  the  partial  derivatives  of  this  expression  with  respect  to  (j)  and  N 
and  setting  the  resulting  equations  equal  to  zero,  the  least  squares  estimates 
are  found  to  be  the  solutions  to  the  following  pair  of  equations: 
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Estimates  -  Least  Squares 

1 


♦is  = 


n 

£  - 

i=l  [N  -  i  +  !]■ 

n  X. 

S  FTTi 

x=i 


(4.58) 


and 


n 


X. 

x 


i-l  (N1S  -  1  +  D: 


n 


X. 


i=l  (NLS  -  i  +  l)2 


n 


(4.59) 


Vi=l  (Nls  -  i  +  iy 


ti=i  (nls  -  i  +  1)*> 


with  the  resulting  estimate  of  the  MTBF  again  being: 

A 

MTBF  (of  the  (j  +  1)  error  occurrence}  =  ■*- 


Z(t.) 


—  -nr 


^LS^LS* 


.  (4.60) 


Tal's32  report  also  provides  estimates  for  <j>,  N,  and  MTBF  based  upon  a  Least 
Squares  Approach  using  the  times  of  error  occurrences,  t^'s,  rather  than  the  time 

between  error  occurrences,  X^s.  It  states  that  the  t's  are  integrals  of  the  X's, 
i  x 

,  and  hence  the  estimates  tend  to  behave  better  as  the  t.'s  fluctuate 


t.  =  T  X. 
i  j 

j=l 


less  due  to  the  cumulative  summary  effect. 

The  estimates  for  <);  and  N  are  derived  by  minimizing  the  sum  of  squared  devi¬ 
ations  of: 


t  (tj  -  expected) *  =  £  (t.  -  £  -Xj  +  TV 


(4.61) 


Taking  partials  and  setting  the  resulting  equations  equal  to  zero,  the  estimates, 
based  upon  the  times  of  error  occurrences,  are  found  to  be  the  solutions  of  the 
following  equations: 


*t  = 


f  A.2 

L  1 

i=l 

V'  t  .A 
Is  i  i 

i=l 


(4.62) 
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with 


and 


with 


1 

A  ’  -  1 

Nt  -  j  +  1 


i=l 


t.B. 
i  i 


i 


E 

j=i 


i _ 

j  +  D2 


(4.63) 


(4.64) 


(4.65) 


Again  the  estimated  MTBF,  after  observing  j  failures,  is: 

A  1 
MTBF  {(j  +  l)  error  occurrence}  =  •* — * -  . 

Vt  - 


(4.66) 


Using  the  various  estimates  of  MTBF,  the  estimated  time  to  remove  the  next  m 
errors,  after  observing  n  failures,  can  be  derived.  Using  any  of  the  previous 
estimates  for  <|>  and  N,  the  estimate  is  obtained  as: 


Estimated  Time  to  Remove  the  Next  m  Errors 


n+m 

53  MTBF  {j  error  occurrence}  (4.67) 

j=n+l 


n+m  . 

53  — — 1 -  •  (4.68) 

j=n+l  <)>(N  -  j  +  1) 


The  only  data  required  for  the  calculation  of  the  estimates  are: 
Data  Requirements 

(a)  Either  the  time  between  error  occurrence  (x^'s)  or 

(b)  The  time  of  error  occurrence  (t/s). 

Once  one  is  recorded,  the  other  (x.  =  t.  -  t.  .)  is  obtained. 

*  11  1“1 
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A  number  of  authors  have  derived  the  large  sample  approximations  to  the  var 
iances  of  these  estimates  (see  References  32  and  33,  also  17,v)  using  the  asympto 
tic  properties  of  the  MLEs.  It  can  be  shown  for  large  n  and  N  that: 


where 


where 


/  A  \  £  _ - _ 

)  =  M  _ (N_-  1  ±. 

j  vMLj  D 

I  s,4 = 


(A. 69 


(4.70 


cov  nml>  ‘•’ml 


n 


(4.71 


n 

E  - 

i=l _ (N 


-  i  +  1): 


j,  -  • 


(4.72 


var  (MTBF^  after  the  nth  error  occurrence}  = 


?  E 


i=l  (N  -  i  +  l)2  (N  -  n)2  (N  -  n)! 


(4.73 


c  .1  r  _ i  _  n  _  zaf 

1  i=l  (N  -  i  +  l)2  (N  -  n)2  (N  -  n)2 


*\ 

■  ►  j-n<|>2(N  -  n)2  + 

j 


+  2<})3(N  -  n)2  (  Y]  x\  +  2E4>3(N  -  n)2  1  -  E24>4 

li=l  ')  i 


(4.7/ 


E  =  £  (n  -  i  +  l)x.  . 

i=l 


(4.7! 


*  There  is  an  error  for  the  variance  of  N  in  this  paper. 
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Using  these  various  estimates  and  variances,  approximate  (1  -  a)  X  100  per¬ 
cent  confidence  intervals  can  be  constructed  for  the  corresponding  population 
parameters  as: 

Confidence  Intervals 


9 

i  V 
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PML 


-  h-a  V 


var!,W'  +«i +  zi-«  Vvat  <w 


KML 


(4.76) 


/  A 

A 

1 

Zl-a 

Vvar*NML}> 

"m.  +  Zl-« 

Vvac  fNML>  / 

(4.77) 

MTBFml  -  Z1-a  Jva7{MTBFML},  MTBF^  +  Z^JvarfMTBF^}  )  *  (4‘78) 

V  2  2  / 

Any  unknowns  in  the  expressions  for  the  variances  are  replaced  by  their  MLEs  and 

a 

Z^_a  is  the  point  taken  from  a  standard  normal  table  such  that  P{Z>Z^_a}  = 

2  ~2 
Large  sample  confidence  intervals  can  be  constructed  using  the  least  squares  esti¬ 
mates.  Schafer,  et.al.  (see  ‘Reference  33)  use  the  result  that  if  Xj,..., 

(4.79) 

X  X 

E{(x.-g.(01,e2))3}  =  P.  (ex,e2)  ,  (4. so) 

which  satisfies  Liapunov's  condition,  that  is, 

n  i  1/3 


X^  are  independent  random  variables  with  finite  moments 
EiX.)  =  (8i,ea),  var{X.)  =  0^(9,, 6S); 


ii?i  PiC0i.e2> J 

lim  - 

n-*»  r  n  ] 

£  o?(9i,e2) 

i-1  J 


1/2 


=  0 


(4.81) 


for  each  finite  value  of  the  unknown  parameters  0*  and  02,  is  three  times  con¬ 
tinuously  differentiable  in  a  neighborhood  of  the  true  values  for  0*  and  02,  then 
the  least  squares  estimators  of  0!  and  02,  obtained  by  minimizing: 


8(0!, 02)  =  Z  [X  -  g  (©! ,02)]2  , 

i=l  1 


(4.82) 
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are  (subject  to  certain  restrictions  on  the  partial  derivatives  of  g^)  asymptoti¬ 
cally  bivariate  normal  with  mean  vector  (0i,02)  and  covariance  matrix: 

i  V ">  .  n  /  3g.\2  n  3g.  9g\. 

isS  /i?!  ai  (set)  !?!  °i  <6l’92)  eel  aet\v- 


E98.  9g  n  /  ogH 

aet  aej  £  “i  (§51 


(4.83) 


where 


.  x  2  2 

n  /8*t\  A  /88i\  /  ",  3Si  S8i 
i=l\36>/  i=l  V89*/  "\i=l  89  ‘  a9* 


(4.84) 
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&  90!  302  \361/ 


(See  Reference  33  for  details.) 

If  the  X.'s  are  the  times  between  error  occurrences,  t 
x 

E*V  =  m  -  (1  -~ot  =  *i«-« 


(4.85) 


(4.86) 


Var{X.}  =  o|  (N,<|>)  = 


<|>2 [N  -  (i  -  l)]2 


with  0i  =  N  and  02  =  4>  in  the  previous  statement  of  the  theorem. 


Thus , 


(4.87) 


<t>[N  -  (i  -  l)]2 


(4.88) 
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and 


*h  Bh 

8<j>  802 


_ -J. _ 

^[N  *  (i  -  1)] 


(4.89) 


If  the  previous  expressions  are  substituted  into  the  asymptotic  covariance  matrix, 
(1  -  a)  X  100  percent  confidence  intervals  can  be  formulated  for  <J>  and  N,  replac¬ 
ing  any  unknowns  by  their  least  squares  estimates.  A  similar  approach  can  be 
taken  to  derive  confidence  intervals  for  4>  and  N  based  upon  the  actual  observed 
times  of  failure,  tJs,  but  it  is  not  developed  here. 

In  seeking  the  estimates  and  their  resulting  confidence  intervals,  the  big¬ 
gest  problem  faced  is  the  difficulty  in  convergence  of  the  numerical  techniques 
employed  to  find  the  MLEs  or  least  squares  estimates.  Difficulties  encountered 
include  (see  Reference  33)  lack  of  convergence,  sensitivity  of  the  iteration 
scheme  to  the  starting  value,  convergence  to  a  saddlepoint  or  invalid  estimate, 
and  nonuniqueness  of  the  estimates.  The  choice  of  a  starting  point  was  especially 
critical  to  the  maximum  likelihood  procedures.  Littlewood  and  Verrall34  have 
shown  that  a  unique  maximum  at  finite  N  and  nonzero  <j>  is  attained  if  and  only  if 


n 

(i  "  1)X. 
- i  > 

(i  >  1) 


otherwise,  the  MLE  for  N  is  “.  Essentially  this  condition  means  that  the  model 
can  only  be  applied  to  software  that  exhibits  software  growth,  i.e.,  X.  >  X.  *. 
In  any  computer  implementation  of  this  model,  the  previous  condition  shoTild  first 
be  verified  to  ensure  a  unique  finite  maximum  exists.  Another  problem  with  the 
MLE  of  N  was  pointed  out  by  Forman  and  Singpurwalla35  concerning  the  instability 
of  the  estimate.  If 


n 


X, 


i=l 


n 


(4.90) 


(i  -  1)X. 


(4.91) 


is  small,  there  is  the  problem  pointed  out  by  Littlewood  and  Verrall,24  but  if  it 
is  large  (so  that  the  times  between  failures  during  the  latter  stages  of  testing 


are  greater  than  the  ones  during  the  earlier  stages),  a  new  problem  with  the  esti¬ 
mate  arises.  The  MLE,  N^,  tends  to  be  close  to  n,  the  number  of  errors  found  to 

date.  This  tends  to  give  a  more  optimistic  view  of  current  reliability  of  the 
program.  It  gives  the  impression  that  the  program  is  very  close  to  being  error¬ 
less  when  in  fact  the  real  error  count  N  may  be  much  larger  than  n. 

To  overcome  this  drawback,  Forman  and  Singpurwalla  suggest  that  the  behavior 
of  the  likelihood  function  be  examined  in  greater  detail.  In  particular,  they 
suggest  the  following  procedure  be  employed  as  a  stopping  rule: 

(a)  Calculate  the  MLE  of  N. 

A  A 

(b)  If  Nml  ~  n,  go  to  step  (c);  if  N^j.  »  n,  observe  another  failure  inter¬ 
val  Xn+1  and  go  back  to  step  (a). 

(c)  Compute 

A 

8(H)  =  (4-92) 

ItNML’ V 

where  L  is  the  likelihood  function 

L(N, <j>)  =  X  $IN  -  i  +  1]  exp  (-<j»[H  -  i  +.l)x.)  (4.93) 

i=l 

and 

<h(N)  =  2 - (4.94) 

1  n  n 

(N  2  x.  -  XI  (i  -  l)x. ) 
i=l  1  i=l 

(d)  Compute 

Formal  (tl)  =  >-«%,  '  N)2/var  ^ml51  (4'95) 

A 

using  the  formula  for  var(N^)  given  earlier. 

(e)  Compare  R(N)  and  RjjqRMAL^  for  various  values  of  N*  If  they  agree 
well,  then  ~  n  is  a  good  estimator  of  N.  If  they  do  not  agree,  then  is 
a  misleading  estimator  of  N;  observe  another  failure  time  interval  and  go  to 
step  (a). 
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These  suggestions  should  also  be  employed  in  any  analysis  utilizing  this 
model. 

In  Forman  and  Singpurwalla1 s  paper,  a  procedure  is  described  for  testing  the 
hypothesis  that  no  errors  remain  in  the  program  and  an  optimal  time  interval  for 
testing  is  developed  based  upon  cost.  The  test  of  the  null  hypothesis 

Hq:  N  =  n,  i.e.,  no  errors  remain  in  the  program 


versus 


H  :  N  >  n 
a 


is  performed  by  exercising  the  software  under  operational  conditi  s  for  an  addi¬ 
tional  b^engtk  of  time.  If  a  failure  is  observed,  the  null  hypothesis  is  re¬ 
jected;  otherwise,  it  is  accepted.  The  additional  length  of  time  is  estimated  as: 


,  -Mi  -  P) 

iMgth  c 


(4.96) 


where  p  is  the  desired  power  of  the  test,  i.e.,  the  probability  of  rejecting  the 
null  hypothesis  when  it  is  false.  The  actually  achieved  power  is  at  least  p 

A 

if  (j)^  is  close  to  (j>. 

The  optimal  aaditional  time,  length*  testing  based  upon  cost  and  mission 

time,  t  ,  is  constructed  as  follows.  If  the  software  fails  during  the  additional 
m 

testing  time,  the  cost  incurred  is 
C i (t )  where  0  <  t  <  t 


length 


(4.97) 


and  Ci(f)  is  a  convex  nondecreasing  function  of  time  representing  the  cost  in¬ 
curred  in  testing  for  time  t.  If  the  software  passes  the  additional  testing  time, 
but  fails  in  operational  use,  the  cost  is 

^'‘length5  *  C*  for  ‘length  i  1  i  ‘length  +  ‘m  -  ”here  (4'98) 

Cg  is  a  fixed  cost  due  to  an  operational  failure  of  the  software. 

If  no  software  error  is  encountered  during  either  testing  or  mission  time,  the 
total  cost  is 


Ci(t,  ..)  for  t,  +  t  <  t. 

1  length  length  m 


(4.99) 


8 

V, 
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The  total  expected  cost  is  therefore: 


/'length  r  'length  m 

MtHe^dt  +  J  (Ci(t-  ri  .,)  +  C2)  <|»e" 


t,  ,  +  t 
length  m 


'length 


length/ 


+  /  c‘tti«*th)*e'*tdt 


t.  +  t 

length  m 


It  can  be  shown  that  when: 


dCi(t) 


>  (J>  C2  exp  {<t»t  }  , 


(4.100) 


(4.101) 


tlength  =  0  min:*-mizes  ttie  expected  cost.  This  means  that  if  the  additional  cost 

of  testing  is  more  than  the  cost  of  an  operational  failure,  no  additional  testing 
should  be  done.  If,  however, 


dCA(t) 


<  4>  C2  exp  {<t>tffl}  , 


(4.102) 


'length 


minimizes  E{C} 


dCi(t) 


t  =  <|>  C2  exp  (<j»tm) 


(4.103) 


A  final  point  concerning  the  Jelinski-Moranda  Model31  is  that  it  is  equiva¬ 
lent  to  Shooman's  Model.  In  Paragraph  4.2.2,  Shooman's  Model  was  given  as: 


(4.104) 


Z(t)  =  K  j  -  ec(T) 


where  a  is  a  proportionality  constant,  1^  is  the  total  number  of  instructions, 
is  the  total  number  of  errors  initially  in  the  program,  and  £c(x)  is  the  cumula¬ 
tive  number  of  errors  corrected  in  the  interval  0  to  t,  normalized  by  the  number 
of  machine  instructions.  Noting  that  £^,  =  N  in  the  Jelinski-Moranda  Model  and 

£  (x)  =  ■— y — ~  for  all  t.  -  <  X  <  t.  ,  (4.105) 

c  I„  l-l  -  -  1 
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the  hazard  rate  function  for  the  Shooman  Model  is  derived. 

_  i  -  1 
‘T 


It  is 


Z(t)  =  K 


N_ 

I. 


I, 


for  t.  ,  <  t  <  t. 
l-l  -  -  i 


*  f-  {N  -  (i  -  1)} 


(4.106) 

(4.107) 


=  4>  {N  -  (i  -  1)1  ,  (4.108) 

letting  <|>  =  y~.  This  is  precisely  the  hazard  rate  for  the  Jelinski-Moranda  Model 
aT  just  considered. 

The  next  few  paragraphs  of  this  report  present,  in  various  amounts  of  detail, 
extensions  that  have  been  made  to  the  basic  model. 


4.2.3. 1  Jelinski  and  Moranda's  Model  1  and  Model  2.  Jelinski  and  Moranda's 
basic  model  cannot  be  applied  to  software  programs  which  are  not  complete.  The 
program  must  be  relatively  stable  with  a  total  of  N  errors  present  initially  in 
the  code.  Their  first  extension  of  the  basic  model36  is  for  programs  that  are 

undergoing  development.  If  at  any  point  in  time  an  error  is  discovered,  an  esti¬ 
mate  of  the  reliability  based  upon  the  percentage  completed  for  the  module  or 
program  can  be  given.  Specifically,  they  let  S(t)  be  the  nondecreasing  fraction 
of  the  total  number  of  statements  which  a  complete  program  has,  measured  at  time 
t,  where  t  is  either  elapsed  wall  clock  or  CPU  time.  Thus  S (0)  =0  and  SO^p) 


=  1, 


where  T^^  is  the  end  time  of  the  program  development. 


The  only  requirement 


about  the  nature  of  the  function  S(t)  is  that  it  be  nondecreasing,  its  values  be 
known  at  the  times  of  error  occurrence  (tj,  t2,...,tn),  and  that  it  be  constant 

during  the  times  between  error  occurrence.  The  hazard  rate  is  then  formulated  as 


Z.(t)  =  (jiS-.j  [N  -  i  +  1]  for  t. 


•1  i  1  i  ‘i 


i  =  1, . . .  ,n 


(4.109) 


Here  N  is  interpreted  as  the  error  content  at  the  end  of  program  development, 
i.e.,  when  S(Tbmt,)  =  1.  S.  .  is  the  fraction  of  the  program  which  was  completed 

cjNL) 

prior  to  the  start  of  the  ith  interval,  i.e.,  =  S(t^_j).  The  likelihood 

function  is  then: 


L (Xj,  X2 


'x-)  -l 


(j)Si_1  [N  -  i  +  1]  exp  {"4»Si_1[N  -  i  +  1]X.}  (4.110) 


where  again 


X. 


x 


t 

i 


t. 

l- 


1 


(4.111) 


The  MLEs  are  obtained  as  the  solutions  to  the  following  system  of  equations. 
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Model  1 

Estimates  -  Maximum  Likelihood 


»ML,1  -  1  +  1 


(4.112) 


E  h-ih 

i=l  1 


E,  S1-1  Vi  ■  1  +  1)Xi  ”  / - 

1~1  „  / 


(4.113) 


i=1\^ML,i  “  x  +  y 

Notice  that  if  S1_1  =  1  for  all  i,  then  the  MLEs  ar  the  same  as  in  the  previous 

paragraph.  The  MLEs  are  very  similar  to  the  ones  obtained  for  the  basic  model  in 
Paragraph  4.2.3.  In  these  equations  there  is  S^X.^,  while  in  the  previous 

section  there  were  X.’s.  The  S.^  has  the  effect  of  reducing  the  time  between 

the  ith  and  (i  -  l)at  error  occurrence  by  the  same  fraction  as  the  percentage  of 
the  program  completed. 

The  MLE  of  the  MTBF,  after  n  errors  have  been  observed,  is  easily  established 


MTBF i  ~ 


(4.114) 


Sftn)(NML,l  “  n^ML,l 


Model  236  also  allows  for  a  developing  program,  but  the  requirement  of  know¬ 
ing  the  fraction  completed  at  each  stage  is  eliminated.  For  a  developing  program, 
it  is  hard  to  envision  a  case  where  the  manager  knows  with  certainty  the  size  of 
the  end  program.  Moreover  the  assumption, 


S(t)  =  S..J 


ti-l  i  1  K  ti  * 


(4.115) 


i.e.,  a  constant  function  between  times  of  error  occurrence,  is  unrealistic.  The 
very  nature  of  a  developing  program  dictates  a  continuously  changing  function  of 

time. 
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For  Model  2,  the  assumption  is  made  that  the  error-making  rate  for  the  pro¬ 
gramming  team  (E^)  is  constant  over  the  time  of  program  development.  The  hazard 

rate  at  time  t  is  then  taken  to  be: 


Z(t)  =  <t>[Gi_1  E  -  (i  -  1)], 


for  t.  .  <  t  <  t. 
i-l  -  -  i 


(4.116) 


'where  <J»  is  the  constant  of  proportionality,  and  is  the  number  of  lines  of 

code  developed  by  the  time  of  the  (i  -  l)st  error  occurrence.  Again  the  basic 
assumption  is  that  the  rate  of  error  occurrence  is  proportional  to  the  number  of 
errors  remaining  in  the  code.  G^  E^  is  the  total  error  count  present  in  the 

G^_j  lines  of  code;  of  which  i-l  have  been  found.  The  likelihood  function  is 

again  expressed  as  (using  the  model  assumptions  of  the  previous  section): 


L(Xl,...,X  )  =  n  4>  [G.  x  E  -  (i  -  1)]  exp  (-p.(G.  ,  E  -(i-l)]).  (4.117) 

i=l  1  p 


The  MLEs  of  (j>  and  E  are  obtained  as  the  solution  to  the  resulting  systems  of 
equations .  " 

Model  2 


Estimates  -  Maximum  Likelihood 


n 

L 


i-l 


2  _  i=l  G,  ,  I  •  (1  •  1) 

♦ml, 2  -  - ii-E - 


(4.118) 


n 


£  Gi-ixi 
i=l 


and 


n 


n  £  G.  .  X, 

n  a  i-l  i 

(G-.i  E  -  (i  -  1))X,  -  1“1 


i-l  p 


u  G 


(4.119) 


i-l 


i-l  Gi_1  E  -  (i  -  1) 


An  estimate  of  the  MTBF  is  then  obtained  as: 

1 


MTBF  2  =  ~- 


(4.120) 


♦mi, 2  tGnEp  -  n) 
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In  the  early  stages  of  development  of  the  program,  the  assumption  that  E^  is 
constant  is  questionable.  As  the  programmers  experience  a  learning  curve  phenom¬ 


ena  the  error  rate  is  expected  to  go  down.  Moreover,  if  the  programmers'  team 
experiences  a  turnover  in  personnel,  with  inexperienced  people  being  hired  at 
various  points  in  the  program  development,  it  seems  hard  to  justify  a  constant  E  . 


The  model  seems  suitable  if  the  development  time  frame  can  be  broken  up  into 
smaller  time  regions  over  which  can  be  taken  to  be  constant. 


4. 2. 3. 2  Lipow's  Extension  Model.  Lipow37  proposed  an  extension  to  the 
Jelinski-Moranda  Model  by  allowing  more  than  one  error  occurrence  during  an  inter¬ 
val  of  testing  and  also  allowing  that  all  errors  found  in  a  given  testing  interval 
need  not  be  corrected  by  the  start  of  the  next  testing  period.  Specifically,  the 
model  assumptions  are: 

Model  Assumptions 

(a)  The  rate  of  error  detection  is  proportional  to  the  current  error  content 
of  a  program. 

(b)  All  errors  are  equally  likely  to  occur  and  are  independent  of  each 

other. 

• 

(c)  Each  error  is  of  the  same  order  of  severity  as  any  other  error. 

(d)  The  error  rate  remains  constant  over  the  testing  interval. 

(e)  The  software  operates  in  a  similar  manner  as  the  anticipated  operational 
usage. 

(f)  During  a  testing  interval  i,  f.  errors  are  discovered  but  only  n.  are 

corrected  in  the  time  frame.  1  1 


The  previous  assumptions  are  identical  to  the  assumptions  of  the  Jelinski- 
Moranda  Model  except  for  (f).  Suppose  there  are  M  periods  of  testing  in  which 
testing  interval  i  is  of  length  x^.  During  this  time  frame,  f^  errors  are  dis¬ 
covered,  of  which  n^  are  corrected.  Assuming  the  error  rate  remains  constant 

during  each  of  the  M  testing  periods  [assumption  (d)],  the  hazard  rate  during  the 
ith  testing  period  is: 

Z(t)  =  <t>[N  -  F±-1]  timl  <  t  <  t.  ;  (4.121) 


where 


<j>  is  the  proportionality  constant, 

N  is  again  the  total  number  of  errors  initially  present  in  the  program. 
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is  the  total  number  of  errors  corrected  up  through  the  (i-l)st 


testing  intervals,  and  t^  is  the  time  measured  in  either  CPU  or  wall  clock  time 


of  the  end  of  the  ith  testing  interval  (x^  =  t^ 


<W- 


The  t.'s  are  fixed 
i 


and  thus,  are  not  random  as  in  the  Jelinski-Moranda  Model.  Taking  the  number  of 
failures,  f^,  in  the  ith  interval  to  be  a  Poisson  random  variable  with  mean 
Z(t.)x.,  the  likelihood  is: 


f . 


M  ((J)  (N  -  F.  , ]x.  ]  1  exp  {-4>[N  -  F.  ,]x.  } 

L(fx,...,fM)  =  n  - -  —  11  1 

i=l 


f . ! 

l 


(4.122) 


Taking  the  partial  derivatives  with  respect  to  4>  and  N  of  £n  L  and  setting  the  re¬ 
sulting  equations  equal  to  zero,  the  MLEs  can  be  obtained  as  solutions  to  the  fol¬ 
lowing  system  of  equations: 


Estimates  -  Maximum  Likelihood 


4>  = 


VA 


NL+  1 


B/A 


and 


M 


M 

-  =  £ 
NL  +  1  -  B/A  i=l 


nl  -  Fi-1 


(4.123) 


(4.124) 


where 


M 

Fm  =  52  f  ,• »  t*ie  total  number  of  errors  found  in  the  M  periods  of  testing, 

i=l  1 


M 


B  =  E  ffi.j  *  «»* 


and 


i=l 


M 

A=  E 

i-1 


x. 

X 


(4.125) 


(4.126) 


the  total  length  of  the  testing  periods.  From  these  MLEs,  the  maximum  likelihood 
of  the  mean  time  until  the  next  failure  is: 


MTBFL  = 


Vnl  -  V 


(4.127) 
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Note  that  if 

f .  =  n.  j=l, . . . ,  M  ;  (4.135) 

J  J  * 

i.e.,  the  number  of  errors  corrected  in  the  ith  interval  is  the  same  as  the  number 
discovered,  then  the  previous  model  reduces  to  a  model  considered  by  Sukert17, 
Gephart  et.al.18  and  Lipow,38 


4. 2. 3. 3  Rushforth,  Staffanson,  and  Crawford* s  Model.  The  last  model  con¬ 
sidered  as  an  extension  of  the  Jelinski-Moranda  Model  is  a  model  by  Rushforth, 
Staffanson,  and  Crawford.39  This  model  was  originally  given  as  an  extension  to 
an  error  generation  model  proposed  by  Shooman.28  A  model  proposed  by  Tal  and 
Barber40  is  also  very  similar  to  the  one  discussed.  The  basic  idea  for  this 
class  of  models  is  to  relax  the  assumption  that  the  error  correction  process  is 
perfect.  This  class  allows  for  the  introduction  of  new  errors  into  the  program 
in  the  correction  of  inherent  ones.  The  specific  assumptions  for  these  models 
are  given  in  the  following. 


Model  Assumptions 


(a)  The  rate  of  error  detection  is  proportional  to  the  current  error  content 
of  the  program. 


(b)  All  errors  are  equally  likely  to  occur  and  are  independent  of  each 
other. 


(c)  Each  error  is  of  the  same  order  of  severity  as  any  other  error. 

(d)  The  error  rate  remains  constant  over  the  testing  interval  of  the  parti¬ 
cular  program  version  undergoing  testing. 

(e)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(f)  No  attempt  is  made  to  correct  detected  errors  at  the  time  of  error 

occurrence.  Instead,  at  specified  points  in  time  t2, . . . ,t. , . . .  a  new  cor¬ 
rected  version  of  the  program  is  provided.  J 

(g)  Of  the  detected  errors  reported,  some  are  corrected,  some  are  not,  and 
some  in  the  correction  process  cause  the  introduction  of  new  errors. 

(h)  The  error  correction  rate,  rc(t),  is  proportional  to  both  the  error 
detection  rate,  rd(t),  and  the  error  backlog  n^Ct),  defined  as  the  difference 
between  the  number  of  errors  detected  by  time  t  minus  the  number  of  errors  cor¬ 
rected  at  that  time.  Specifically,  rc(t)  is  taken  as: 

rc(t)  =  a  rd(t)  +  pnb(t).  '(4.136) 
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(i)  The  rate  of  error  generation  is  assumed  proportional  to  the  error- 
correction  rate. 

r£(t)  =  Y  rc(t).  (4.137) 

(j)  The  error  detection  process,  nd(t),  is  completely  known. 

From  these  assumptions,  Rushforth,  Staffanson,  and  Crawford's  Model  is  formu¬ 
lated  as  shown  in  the  following. 

For  any  model  version  j , 

Nj  =  N  -  nc(tj_j)  +  n£(tj_p  ,  (4.138) 


where 


N.  is  the  total  number  of  errors  present  in  version  j, 
J 

N  is  the  initial  number  of  errors  present, 


n  (t.  „)  is  the  total  number  of  errors  corrected  up  through  program  version 
cv  j-17 

j-1,  and 

n  (t.  ,)  is  the  total  number  of  errors  introduced  into  the  program  in  the 

£  N  J  -1 

correction  of  the  previous  (j  -  1)  versions.  From  assumption  (a), 


r„(t)  =  0  N. 


for  all  t's  ti_1  <  t  <  t± 


(4.139) 


where  4*  is  the  proportionality  constant.  From  assumption  (h),  if  nc(t)  is  the 
number  of  errors  detected  up  through  time  t,  then, 

n^t)  =  nd(t)  -  nc(t)  (4.140) 

so  that: 

rc(t)  =  a  rd(t)  +  p  (nd(t)  -  nc(t)).  (4.141) 

Finally,  using  the  fact  that: 

n  (t)  =  f  r  (x)dx  (4.142) 

e  o  e 

and 

nc(t)  =  f  rc(x)dx  ,  (4.143) 
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it  can  be  drawn  from  assumption  (i)  that: 


n  (t)  =  f  r  (x)dx  =  Y  f  rc(x)dx  =  Y  nc(t). 

0  0 

(4.144) 

Hence, 

rd(t)  =  <J>N. 

(4.145) 

=  4>  [N  -  +  ne  (tj.j)) 

(4. 146) 

=  (J)  [N  -  n^t-.p  +  V  nc  (tj_1)) 

(4.147) 

=  4>  (N  -  (1  -  Y)  nc  (t.^))  tj.j  <t<t. 

(4,148) 

and 

rc(t)  =  a  rd  (tj.,)  t  -  ^(tj.j)) 

(4.149) 

=  a  (j)  Nj  +  0  [nd  (tj-;l)  -  nc  CtJ_1) ] 

=  a  <t>  [N  -  (1-Y)  nc  (tj_1) ]  +  p  [nd  (t^)  -nc  ] 

(4.150) 

i  1  i  t j 

(4.151) 

If 

4>a  =  (1  -  YH 

(4.152) 

and 

Na  =  N/(l  -  V)  , 

(4.153) 

the  two 

equations  become: 

rd(t)  =  *a  1Na  '  nc 

(4. 154) 

and 

rc(t)  =  «  ♦,  [Ha  -  nc  (tj.j)]  +  P  Ind  (tj.i)  '  "c  (t..,)] 

(4.155) 

involving  five  unknowns  N  ,  <j>a,  a,  P,  rc(t).  From  assumption  (j)>  n^Ct)  is  known 
exactly  so  that, 


rd(t)  = 


dnd(t) 

dt 


(4.156) 


is  known  exactly. 


Rushforth,  Staffanson,  and  Crawford  show  through  a  linear  system's  approach 
that  the  two  equations  can  be  expressed  as: 

n,(t.)  =  -  *a  +  Na  »a  At.  (4'157) 


and 


ac(t,)  =  +  nc(t..j)  -  nc(t^_j)  A  IP  +  «♦,! 

j  J 


with 


+  “Na  ♦a  At. 


\  -  tj  -  tj.,  • 


(4.158) 


(4.159) 


Here  there  are  two  equations  in  four  unknowns  Nfl,  $a,  a,  and  p.  If  the  addi¬ 
tional  assumption  is  made  that  A^.  is  a  constant  T  (i.e.,  the  new  program  versions 

are  introduced  as  equally  spaced  points  in  time),  then  the  two  equations  can  be 
reexpressed  as: 


where 


and 


^nd(t.)  -  nd(t._j) 


6(t.)  =  |  "  y  3‘ 

lnc(V  -  “c 


=  v.  Li'lB 


2x1 


C  ■*'  ) 

\pT  1  -  TIP  +  <t>a0f]y 


(4.160) 


(4.161) 


(4.162) 


Their  paper  develops  the  estimates  of  the  four  unknowns  using  a  least  squares 
approach.  The  quantity  for  which  the  minimum  is  sought  is: 

K 

£  (6(t.)  -  observed  increments)2  (4.163) 

i=l  ~  1 


where  X  is  the  number  of  program  versions, 


N  <t> 


(4.164) 
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and  the  observed  increments  = 


M  u  * 

r;  Csj 


'W  •  “d^j-iV 

^nc(tj)  - 

difference  in  number  detected  from  one  version 
'to  the  next 

difference  in  number  corrected  from  one  version 
to  the  next 


(4.165) 


/ 


A  computer  program  was  developed  to  perform  this  nonlinear  minimization  process 
and  is  provided  in  their  paper. 


4.2.4  Schick-Wolverton  Model 


;  a 


The  next  class  of  model  considered  was  originally  proposed  by  George  Schick 
and  Ray  Wolverton.6  Their  model  assumes  that  the  hazard  rate  function  is  not 
only  proportional  to  the  number  of  errors  in  the  program,  but  proportional  to  the 
amount  of  testing  time  as  well.  Their  logic  is  that  as  testing  progresses  on  a 
program,  the  chance  of  detecting  errors  increases  because  of  ”zeroing-in"  on 
those  sections  of  code  in  which  errors  lie.  Specifically,  their  model  is  based 
on  the  following  assumptions: 

Model  Assumptions 

(a)  The  rate  of  error  detection  is  proportional  to  the  current  error  content 
of  a  program  and  to  the  amount  of  time  spent  in  testing. 


other. 


(b)  All  errors  are  equally  likely  to  occur  and  are  independent  of  each 


(c)  Each  error  is  of  the  same  order  of  severity  as  any  other  error. 

(d)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(e)  The  errors  are  corrected  instantaneously  without  introduction  of  new 
errors  into  the  program. 

The  one  major  difference  between  these  assumptions  and  the  Jelinski-Moranda 
Model  is  assumption  (a)  with  the  error  rate  also  being  proportional  to  the  amount 
of  testing  time. 

The  form  of  the  hazard  rate  function  is: 


Z(X.)  =  4*  [N  —  (i  -  1)]X.  , 


(4.166) 
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where  X^  is  the  amount  of  testing  time  spent  between  the  discovery  of  the 
(i  -  l)st  error  at  time  t^_1  and  the  ith  error  at  time  t^.  The  quantity  <j>  is  the 

proportionality  constant  of  assumption  (a)  and  N  is  the  total  number  of  errors 
initially  in  the  program.  Figure  4-2  is  a  plot  of  this  function  over  time. 
Using  the  relationship  established  earlier  between  the  hazard  rate  function,  the 
reliability  function,  and  the  MTBF,  it  can  be  s*en  that: 

2 

X. 

R(X^)  =  exp  {-<)>  [N  -  (i  -  1)]  2^  }  (the  Rayleigh  distribution)  (4.167) 


t-w, 


and 

MTBF=  /  ROq)^  = 


(4.168) 


To  develop  the  MLEs  of  <}>  and  N,  suppose  errors  are  discovered  at  times  ti, 
. *.,tn  and  suppose  X.  =  t^  -  t^_^.  Then,  from  equation  4.167. 

.  2 


RCX.^)  =  exp  { -<J) [N 


X. 


(i  -  1)1  ! 


(4.169) 


Now  f(X^)  =  -R'(X^)  so  the  distribution  of  the  time  between  the  (i  -  l)st  and  ith 
error  is: 

2 

X. 

f(X.)  =  4>[N  -  (i  -  1)]X,  exp  H[N  -  (i  -  1))  ~  }  .  (4.170) 

z(t) 


FIGURE  4-2.  SCHICK-WOLVERTON  HAZARD  RATE  FUNCTION 
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Using  the  asymptotic  properties  of  the  MLEs,  it  can  be  shown  (References  17 
and  33)  that  for  large  n,  the  large  sample  variances  are: 


S3 

5  *i 

i 

v  ’ 

I 

i 

*,N 

i  * 

* 


var  *<W  s 


var  ^Nsw^  = 


n  1 

r  — 1 - 

i=l  (Ngw  -  i  +  l)2 
_ 


4>2d3 


cov  (Nsw,  (j.sw)  = 


n  2 
2DX 


and 


n 


i=l  (Ngw  -  i  +  l)2 


(4.177) 


(4.178) 


(4.179) 


(4.180) 


100  X  (1  -  a)  percent  confidence  intervals  can  then  be  constructed  a 3  before  as: 
Confidence  Intervals 


*sw  ■  zi-«  «  •  »SW  +  Zl-oVV3r  (W 


and 


where 


Jl-a 

2 


NSW  ‘  Zl-a  Vvar *NSW*  ’  NSW  +  Zl-a  Vvar 


is  chosen  from  a  standard  normal  table  so  that: 

F  (Z  >  Z^}  =  £  . 

2 


(4.181) 


(4.182) 


(4.183) 


The  expected  time  to  remove  the  remaining  N  -  n  errors  after  n  errors  are 
discovered  is: 


Expected  Time  to  Remove  the  Remaining  N  -  n  Errors 
N 


■  £  J 

j=n+l  1 


n 


2<j»tN  -  (  i  -  1)  ] 


(4.184) 


(3 


>/>,  *1 
M 
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(4.185) 


The  large  sample  MLEs  of  this  parameter  are  therefore: 
Estimated  Time  to  Remove  the  Remaining  N  -  n  Errors 


_  Nsvfn 


2<j)sw1 


(4.186) 


As  in  the  Jelinski-Moranda  Model,  the  least  squares  estimates  of  4>  and  N  are 
obtained  by  minimizing: 


£  (X.  -  MTBF . ) 2  =  £ 

i=l  1  1  i=l 


(4.187) 


Taking  the  partial  derivatives  with  respect  to  4>  and  N,  the  least  squares  esti¬ 
mates  are  obtained  as  the  solutions  to  the  following  pair  of  equations: 


Estimates  -  Least  Squares 


*XS,SW 


l  NLS,SW  "  1  +  1 

't  -* . — ,v 

1=1  ^NLS,SW  “  1  +  ^  j 


(4.188) 


£  —*  2/2 

-1  <NLs,sw“i  +  1) 


2<Jils,sw  i=1  (nls.sw  "  1  +  1)2 


.  (4.189) 


The  asymptotic  variances  for  these  estimates  can  be  developed  using  the  approach 
taken  in  the  Jelinski-Moranda  section  of  this  paper. 

For  this  report,  only  one  basic  extension  of  the  Schick-Wolverton  Model  pro¬ 
posed  by  Lipow38  is  considered. 


4.2.4. 1  Lipow1 s  Extension  to  the  Schick-Wolverton  Model.  Lipow 's  Model 
uses  the  same  assumptions  as  the  Schick-Wolverton  Model  except  assumption  (a)  of 
the  previous  paragraph  is  replaced  by 
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6SS  l-  ' 

[V^N 

Cvi 

g  m 


(a)  The  rate  of  error  detection  is  proportional  to  the  current  error  con¬ 
tent  of  the  program  and  the  total  time  previously  spent  in  testing  including  an 
"averaged"  error  search  time  during  the  current  time  interval  of  testing. 

Another  way  in  which  this  model  differs  from  the  previous  one  (hence,  it  is 
not  a  true  extension)  is  that  it  is  an  error  count  model  rather  than  one  using  the 
"time  between  error  occurrences."  Suppose  f.  errors  are  discovered  during  the 

i 


ith  testing  interval  and  suppose  F^  =  ^  f ^ 
discovered  up  through  i  testing  intervals. 


is  the  cumulative  number  of  errors 


that: 


Based  upon  the  model  assumptions  for  the  hazard  rate  function,  it  can  be  seen 


Z(x.)  =  4>(N  -  Fi_1  ]  [X..J  +  x.]  , 

T 


(4.190) 


where 


x^  is  the  amount  of  testing  time  spent  between  the  end  of  testing  period 
i  -  1  and  the  end  of  the  ith  testing  period, 


a  «  ''V> 


ja  m 


R  i 

f’i' 


%  m 


ft  m 

r,  I  [■< 

ft 


h-i =  £  xj  > 

i.e.,  the  cumulative  amount  of  testing  time  spent  through  i  *  1  intervals. 
Since  the  reliability  function  is  related  to  the  hazard  rate  as: 

x 

R(x)  =  exp  {-  f  Z(v)dv}  , 


(4.191) 


(4.192) 


* 

A  >>:'• 
A:-' 


1  h 

!-V0 

S  & 

/  *1  • 


R(x.)  =  exp  {  f 1  4>[N  -  Fi_13  [X._1  +  ^]dv}  ; 


R(xt)  =  exp  {~4>[N  -  Fi_1]  [Xi_1x.  +  x?]  } 

r 


(4.193) 


i94) 


a  m 

<  -1  L  W 


'm 


@  m 

L"»J  L.o 

av 
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Also,  as  noted  earlier, 


00 

MTBF  =  /  R(x.)dx. 


(4.195) 


im 


$ 


MTBF 


00  X? 

=  /  exp/^IM  -  Fj.jl  [X^  +  ji  ]}  dx1 
o  l 


0#  f  ,  (N  -  F.  ,] 
r  exp^-f  _ _ i-l 


(4.196) 


[*?  +  «i-i*i  ‘  «?.,  -  «£-:]}<*! 


00  f  (N  -  F.  ,]  ,  7.  2  nl 

=  /  exp|-<ti  - 5-^-  |_(x.  +  2X..1)  -  4Xi_1J|dx. 


(4.197) 


(4.198) 


=  exp  {*  [N  -  F-.j)  X^j}  /expj-*  [N  "JW  (xt  +  2X._1)2Jdxi 


(4.199) 


exp  |<|»  [N  -  F.^IX^ 


00  , 

, _ _  r  _ L 

■i-i]  o  I  5f 


r  lx,  -  <-2x..1)i°'| 

expj-  -5-7 - 5 - T  !■  dxi 


( ?l »  -\-y 


exp  {<»[N  -  F,  ,]X*  }  ^  :  -  •  1/2 

1  1  1  1  ilW  -  Fj.j] 


,  =  exp  i*[S  -  F.  ]X|  ). 

♦  IH  ‘  F-j] 


(4.200) 


(4.201) 


(4.202) 
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Sukert17  gives  the  MLEs  for  <j>  and  N  as  the  solutions  to  the  following  system  of 
equations : 

Estimates  -  Maximum  Likelihood 


£  (i*LSW  '  + 


(4.203) 


M  f.  a  M  x. 

£  ”  — -  -  *I,SW  £  xi  (-i-l  +  T]  ■  (4-204) 

11  klsw  '  Fi-1  1-1 

A 

(Note:  Sukert  left  out  ^j^SW  '*'rl  ^is  e(l'a&^on  corresponding  to  4-204.) 

M 

M  is  the  number  of  testing  intervals  and  FM  =  £f.  t^le  total  number  of  errors 

i=l 

found  in  all  the  testing  intervals.  The  MLE  of  the  MTBF  and  the  expected  time  to 
remove  the  remaining  N  -  F^  errors  are: 


MTBF  = 


^LSW^LSW  “  FM] 


exP  *LSW  lNLSW  '  F1 


V3} 


(4.205) 


Estimated  Time  to  Remove  Remaining  N  -  FM  Errors  =  . 


g1  yn  fpjjw  ^lsw-f^ 

J=Fw  /  A  A 

VWlsw'V 


(4.206) 


Sukert.  also  provides  the  large  sample  asymptotic  variances  of  the  estimates 
for  4>  and  N  as: 


^T.sw^  =(  53 


i=1  (NLSW  “  Fi-1^ 


(4.207) 


var  {nlsw}  =  To  Z 

*ls\P 


(4.208) 


me 
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where 


c°v  <NLSW>W  = 


£ *i  (?!-!  +  !i> 
i=l  1  1  1  2 

D 


^LSW  '1=1  ^NLSW  "  Fi-1)2/ 


£  *i  +  rj  • 


(4.209) 


(4.210) 


As  in  the  previous  paragraphs,  large  sample  confidence  intervals  can  be  con¬ 
structed  for  <j>  and  N.  From  the  previous  formulas,  it  can  be  noticed  that  if: 

X.  .  is  set  =  0  and  for  all  i, 

-l-l  ’ 

f^  =  1  i  =  l,...,ri,  so  that 
Fi-1  =  i  “  1’ 

then  the  formulas  for  the  estimates  and  their  variances  become  those  associated 
with  the  Schick-Wolverton  Model. 


4.2.5  Generalized  Poisson  Model 

The  Generalized  Poisson  Model  (GPM)  was  given  in  a  report  by  Schafer,  Alter, 
Angus,  and  Emoto33  for  the  Hughes  Aircraft  Company  under  contract  to  the  Rome  Air 
Development  Center.  Their  model  can  be  considered  to  be  analogous  in  form  to 
both  the  Jelinski-Moranda  and  Schick-Wolverton  Models  but  taken  within  the  error 
count  framework.  With  a  slight  modification,  it  can  be  shown  to  be  an  extension 
of  Lipow's  Model  as  well.  The  model  assumptions  are  given  in  the  following. 

Model  Assumptions 

(a)  The  expected  number  of  errors  occurring  in  any  time  interval  is  propor¬ 
tional  to  the  error  content  at  the  time  of  testing  and  to  some  function  of  the 
amount  of  time  spent  in  error  testing. 

(b)  All  errors  are  equally  likely  to  occur  and  are  independent  of  each 
other. 

(c)  Each  error  is  of  the  same  order  of  severity  as  any  other  error. 

(d)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(e)  The  errors  are  corrected  at  the  ends  of  the  testing  intervals  without 
introduction  of  new  errors  into  the  program.  (Note:  Errors  discovered  in  one 
testing  interval  can  be  corrected  in  others;  the  only  restriction  is  that  the 
error  corrections  come  at  the  end  of  the  testing  intervals.) 
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Using  the  assumptions,  their  model  is  constructed  as  follows.  Suppose  the 
testing  intervals  are  of  length  xi,...,xn  and  suppose  f^  errors  are  discovered  in 

the  ith  interval.  At  the  end  of  the  ith  interval,  a  total  of  M.  errors  are  cor- 

rected.  In  the  previous  extensions  of  the  Jelinski-Moranda  and  Schick-Wolverton 

i 

Models,  M.  was  set  to  £  f.,  i.e.,  all  errors  found  in  an  interval  are  corrected 
1  j=l  J 

at  the  end  of  the  interval.  This  model  relaxes  this  assumption. 


From  assumption  (a) , 

”  M._1]gi(x1,x2,...,xi) 


(4.211) 


where  <j>  is  the  proportionality  constant,  N  is  the  initial  number  of  errors,  and  g^ 

is  some  function  of  the  amount  of  testing  time  spent  previously  and  currently. 
Usually,  g^  is  nondecreaaing  with  the  logic  that  as  more  time  is  spent  in  testing, 

more  errors  are  discovered.  In  the  paper  by  Schafer,  et.al.,  the  function  g.  is 
taken  to  be:  1 

8i^i,x2,. . .  ,xi)  —  x^  .  (4.212) 

This  restriction  is  relaxed  to  show  a  broader  class  of  adaptability. 

For  example,  if 

8i(xi,...,xi)  =  xi  ,  (4.213) 

then  the  resulting  formulas  for  the  estimates  are  the  same  as  the  Jelinski- 
Moranda  Model;  if 


8i(xi,...,xi)  r  x?/2, 


then  the  formulas  are  the  same  as  the  Schick-Wolverton  Model;  and  if 


.  .  ,x.)  =  x ,  J  £  x.  +  ~ 
1  I  j=l  J 


(4.214) 


(4.215) 


Lipow's  formulas  are  obtained. 

From  assumptions  (a)  and  (b) ,  the  joint  density  of  the  f/s  is: 


f(fx,...,f  )  =  n  f(f  ) 

n  i=i  1 

n  [<j>(N  -  M.  ,)g.(x1,...,x.)]fi 

nr  111  1 


(4.216) 


exp  {-<j>(N  -  Mi_1)gi(x1,...,xi)}; 


(4.217) 
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i.e.,  f^  is  Poisson  with  mean  =  <j>(N  - 

Hence,  the  likelihood  function  L  is: 


(4.218) 


so  that 


L(4>,N)  =  n  f(f.) 

i=l  1 


£nL(<J>,N)  =  £  £n£(f ,) 

i=l  1 


(4.219) 


(4.220) 


tn  n 

f  ,£n(j)  +  Vf£n(N-M.,)  +  £f.£ng.  (4.221) 

i=l  1  i=l  1  11  i=l  x  1 

-<J>  £  (N  -  M.  ,)g.  -  f  £nf .  I  . 
i=l  111  i=1  i 


Thus,  taking  the  partial  derivatives  with  respect  to  <{>  and  N,  the  following  is 
obtained: 


f.  n 


-  -  £  » '  Mi-i>«i 

n  f .  n 

-  S  W^TTJ  -  ♦  .?  % 

1=1  N  l-l7  1=1 


(4.222) 


(4.223) 


Thus  the  MLEs  of  <}>  and  N  are  solutions  to  the  following  pairs  of  equations: 

i 

Estimates  -  Maximum  Likelihood 


n 

£** 


N™  £  *i  -  E 

1=1  i=i 


(4.224) 
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n 

£  Si 

i=l 


'GPM 


n 

E 


i=l  (N, 


GPt 


Xi  U.1  xi 

Note  if  gi  is  the  appropriate  function  "  2“  >  or  xi  ^  xj  +  * 


and  if 


a  n  g  ^  £i 

*GPM  S  (NfiDM  -  H,.,)'  =  i=l  (SGPM  -  ' 


(4.243) 


then  not  only  are  the  MLEs  for  the  GPM  the  same  as  the  previous  paragraphs,  but 

the  asymptotic  variances  are  the  same  as  well. 

,  \ 

From  the  GPM  formulation,  the  expected  number  of  errors  in  the  (n  +  l)st  in¬ 
terval  of  testing  is 

E{  £n+l|  =  ♦  01  '  VSn+l  (X1 . Vl>  '  (4’2W) 

where  x  ,  is  the  anticipated  testing  time  for  the  (n  +  l)st  session.  The  MLE  is 
n+1 

therefore 


E{fn+l}  =  ^GPM  (NGPM  "  Mn)gn+1  (xi  »  *  *  * ’ Vl}  * 
Its  asymptotic  variance  is  therefore 

r ; r  fa 

var{  E{fn+lj)  ~  3* .  9H  6  1  8 


(4.245) 


8EU. 


(Ngpm  -  V8n+1  *GPM  8n+l 


^NGPM  ~  Mnfa+1  ^GPM  8n+l 


(4.246! 

* 

(4.247! 


~  (NGPM  '  M»)8n+1  ^GPM  8n+l 


^GPM  y  gl _ 

^  »  <«OPK 


n 

*  S  g 
i=l  8i 


t  »4 
1=1  1 
“D* 


n 

S'* 


^NGPM  "  Mn^  8n+l 

A 

s  *GPM8n+l 


(4.248 
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In  the  Hughes  report,  the  least  squares  estimates  are  also  derived.  Those 
estimates  are  chosen  to  minimize: 


n 


S  =  £  [f\  -  <j>  (N  -  M.^g.}2. 

i=l 


By  taking  the  partial  derivatives  of  S  with  respec  t-  to  <j>  and  N, 

n 


(4.253) 


and 


as 

8<j) 


as 


*  -2  £  iv  -  ViV1  «  ■ 


(4.254) 


n 


3N  ■  ‘2  £  ‘V  ♦»  - 

1=1 


(4.255) 


The  least  squares  estimates  are  then  obtained  as  the  solutions  to  the  following 
system  of  equations: 

Estimates  -  least  Squares  v  > 

n 


4> 


£ 

i=l  1  1 


LS,GPM 


£  (N  -  M  )g2 
1=1  i-1  i 


(4.256) 


and 


£  Vi.  « 


n 


i=l 


^LS , GPM  “  Mi-P  "  ^LS , GPM  £.  ^NLS,GPM  ’  Mi-l^8i 
>  >  1=1  > 


=  0 


(4.257) 


Using  the  large  sample  results  for  leasft  squares  estimators  given  in  Para¬ 
graph  4.2.3  of  this  report  and  using  the  Hughes  report,  it  can  be  shown  that: 


m 


Sis, GPM  =  Vs, GPM  ’  NLS,GPm' 


(4.258) 


0 
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is  asymptotically  normally  distributed  with  mean  vector  (4>,N)  and  covariance  na- 
trix: 


/i?i  ^S,GPH  <NLS,GPM  ‘  ♦L.GPM^LS.CPM  *  Hi-p2gi^ 


^CPK 


:  J-  £ 

A2  1 


IE 


*LS,GPM(HLS,GPH  *  J-j  ♦lS.GPH  (NLS,GPH  “  Mi-Pgiy 


(4.2159) 


where 


£ 


l  - 


_  A  Jl  A  A 

^  ^LS ,  GPM8i  “  ^LS ,  GPM^NL5 ,  GPM  “  Mi-Pgi 


(4-2<>0) 


n 


fcj 

and 

*v« 

I 

A  = 

I 

R; 

=  (£ 
U*i 


A2 

*LS,GPMgi 


^LS.GPM^LS.GPM  '  Mi-P8i  (NLS,GPM  ‘  Mi-l)2gi> 


^  (NLS,GPM  "  Mi-l)2gi^“  *LS,GPM^NLS,GPM  "  Mi-l)gi^ 


(4.2<>1> 

\ 

Large  sample  confidence  intervals  can  then  be  constructed  using  tbetse  results 


The  data  requirements  necessary  to  implement  this  model  are: 

Data  Requirements 

(a)  The  Lengths  of  the  various  testing  intervals,  i.e.,  xi,...,xn> 

(b)  The  number  of  errors  corrected  at  the  end  of  each  testing  period,  and 
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(c)  The  number  of  errors  discovered  in  each  interval  of  testing,  i.e.,  f.'s. 

In  the  Hughes  report33,  an  extensive  analysis  was  done  on  the  properties  of 

the  least  squares  and  MLEs  for  the  specific  model  of  g^(xj, . . .  ,x^)  =  xi  with  a  a 

positive  integer.  Thus  the  GPM  could  be  taken  as  a  three  parameter  estimation 
problem  of  N,  0,  and  a  or  a  two  parameter  problem  if  a  is  specified.  The  report 
noted  problems  with  lack  of  convergence,  convergence  to  nonoptimal  solutions,  and 
lack  of  uniqueness  for  the  estimates.  The  report  noted  that  the  MLEs  had  a 
greater  tendency  to  these  problems  than  the  least  squares.  Both  procedures  were 
dependent  upon  a  '‘good"  initial  starting  value  for  the  estimates  in  order  to 
achieve  convergence.  Other  problems  with  the  estimation  procedures  included  ob¬ 
taining  solutions  which  violated  model  assumptions  and  oscillation  of  the  esti¬ 
mates  in  the  convergence  process. 


4.2.6  Geometric  Model.  This  model  was  proposed  by  Moranda  (References  41  and  42) 
and  is  a  variation  of  the  Jfelinski-Moranda  "De-Eutrophication"  Model.  It  is  an 
interesting  model  because,  unlike  all  of  the  previous  models  discussed,  it  does 
not  assume  a  fixed  finite  number  of  errors  in  the  program,  nor  does  it  assume  the 
errors  are  equally  likely  to  occur.  It  assumes  that  as  debugging  progresses,  the 
errors  become  harder  to  detect.  By  operating  on  the  premise  that  a  program  is 
never  completely  error  free  (because  of  error  introduction  in  the  process  of  cor¬ 
recting  a  detected  error),  this  model  can  be  utilized  for  error  analysis.  The 
specific  model  assumptions  are: 


Model  Assumptions 

(a)  There  is  an  infinite  number  of  total  errors  (i.e.,  the  program  is  never 
error  free) , 

(b)  All  errors  do  not  have  the  same  chance  of  detection, 


(c)  The  detections  of  errors  are  independent, 

(d)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage,  and 

(e)  The  error  detection  rate  forms  a  geometric  progression  and  is  constant 
between  error  occurrences. 


From  these  assumptions,  the  hazard  rate  for  this  model  is  of  the  form 

Z(t)  p  D^1”1  (4.262) 

for  any  time  t  between  the  occurrence  of  the  (i  -  l)st  error  and  the  ith.  The 
hazard  rate  function  is  initially  a  constant  D  which  decreases  in  a  geometric 
progression  (0  <  <J>  <1)  as  error  detection  occurs.  Figure  4-3  is  a  graphic  repre¬ 
sentation  of  this  hazard  rate  function. 
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FIGURE  4-3.  GEOMETRIC  ”DE~EUTRpPHICATION”  PROCESS 


Notice  from  the  graph  the  ratio  of  the  change  in  the  error  detection  rate, 
2(t) >  from  the  discovery  of  the  ith  error  to  the  (i  +  l)st,  i.e., 


Change  in  Z(t)  on  the  discovery  of  the  ith  error 


i-1  i-2 

D<t>  -  D<t> 


Change  in  Z(t)  on  the  discovery  of  the  (i+l)st  error  nJ  nAi"l 


D^-Df 


=  ?>l 


(4.263) 


Thus  the  size  of  the  step  gets  smaller  as  errors  are  discovered.  This  means  that 
latter  errors  are  more  difficult  to  find  and  do  not  have  as  dramatic  effect  on 
lowering  the  error  rate  as  earlier  detected  ones. 


Again,  if  X^  =  t^-  t^_^  is  the  time  of  discovery  between  the  ith  and  (i  - 


l)st  error,  then  using  assumptions  (c)  and  (e),  the  X/s  are  assumed  independent 


exponentials  with  rate  Z(t^),  i.e., 


f(X.)  ==  D4>i_1exp  {-D^X.}  . 


(4.264) 


The  likelihood  function  for  the  X.'s  is  then: 

i 


n 


n 


L(X!,...,Xn)  =  n  f(X. )  =  Dn  n  <|>i"1exp  {-D  £  ♦1"ix.} 

n  i=l  1  i=l  i=l  1 


n 


i-1. 


(4.265) 
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Thus  the  log  of  the  likelihood  function  is : 


n  n  .  , 

=  a£nD  +  E  (i-l)2n<j>  -  D  T 

i=l  i=l  1 


(4.266) 


The  MLEs  of  D  and  $  are  obtained  as  the  solutions  to  the  following  pair  of  equa¬ 
tions: 


32nL  _  n  v->  j.i“lw  1  _  * 

-W  -  5  £  *  xi  -  0 

1=1 


(4.267) 


8lnL  _  A  (i  -  1) 

3<t>  ”  <♦> 


-  D  E  (i  -  l)^i-2X.  =  0. 

i=l  1 


(4.268) 


Estimates  -  Maximum  Likelihood 


That  is: 


D  = 


E  *cxi 

i=l  b  1 


n 

n  a 

E  <t»GXi 

i=l  1 


n+1 
2  • 


(4.269) 


(4.270) 


From  the  estimates,  the  MLE  of  the  MTBF  after  n  errors  have  been  observed  can 
be  estimated  as: 


MTBF  =  E  X 

*  n+1*  n  An 

DG 


(4.271) 


The  model  cannot  be  used  to  estimate  tht  total  number  of  errors  in  the  program 
but,  it  can  be  used  to  estimate  the  "purity"  level  after  n  errors  are  observed. 
The  estimated  degree  of  "purification"  for  a  program  is  usually  given  by  the 
ratio: 


Z(t  )  -  Z(t  )  n  n.n 

—2 - ai  =  e^L  «  ! .  „■>  , 

2(t0) 


{i*. ill) 
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the  change  in  the  hazard  rate  function  from  the  beginning  of  testing  to  the  end 
versus  what  it  was  at  the  beginning.  The  estimate  of  this  purification  level  is 
therefore:  ... 


n 


PL  =  1  - 


(4,273) 


The  large  sample  estimates  for  variances  of  these  estimators  are  derived 
again  using  the  large  sample  properties  "of  the  MLEs  (see  Paragraph  4.2.5.). 


Now 


and 


Hence , 


and 


3  2dL  _  _  n_ 
3D2  D2 


2 

3  £nL 
3<t>3D 


32AnL 


“  £  (i  " 

i=l 


(i-2) 


X. 

l 


n-2 


«  -  n<n  ~  U  -  d  £  t(i  +  i)  4ii_1x. 

_  .  n  1 


2(j»i 


i=l 


i+2 
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1  D2 

(4.277) 
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(■ij 

f  32Aul) 
[  3<t»3D  J 

►  =  t  (i  -  D<r2E 

i=l 

txi> 

(4.278) 

i 

=  £  U  -  1H1’2 

i=l 

t1'1) 

(4.279) 

■if  £  Ci- 

_  n(n  -  1) 

“  2D(j» 

(4.280) 
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L*..< 

f  3  AnL 
1  3<j»2  . 

L  n(n  -  1)  +  D  ‘jj* 

)  2(j)2  i=l 

i(i  +  1H1’1  E  {xi+2} 

(4.281) 
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n(n  -  1)  +  1  T(n_ 

202  02L 

g(2n  -  l)(n  -  1) 


2) (n  -  l)(2n  -  3)  .  (n  -  2) (a 
S  +  2 


-  1) 


51 


(4.283) 

(4.284) 


Hence  for  large  n, 


*6  = 


(4.285) 


is  asymptotically  bivariate  normal  with  mean  vector  (D,0)'  and  covariance  matrix 

i-1 


E 


n(n  -  1) 
2D0 


2D2(2n  -  1) 

n(n  +  T) 

-6D0 
n(n  +1) 


n(n  -  1) 

2D0 

n(2n  -  l)(n  -  1) 

602 

-6D0 
n(n  +  1) 


(4.286) 


n(n 


1202 
-  l)(n  + 


(4.287) 


i.e. 


and 


a  Mr(2a  ~  1) 

variDoi  ~  Hfornr 


a  ^  ,1202 

var ~  iy(a  +  l) 


A  A  6Dr,0« 

Cov(W  s  'nln-hr 


(4.288) 


(4.289) 


(4.290) 
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The  estimated  variances,  of  the  MTBF  and  PL  are  then  obtained  from: 


var 


s«»F)  =(|g  ,  if)  (|f  ,  If) 
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where 


8  =  - 


(4.292) 


where 


var  {PL)  ~  (§§ 


f  =  1  -  <j>“  . 


These  can  be  shown  to  be: 


i)  £.  0 


(4.293) 


(4.294) 


var(MTBF)  =  ,.  -  -1) 

Dg  <|>g  n(n  +  1) 


var{PL}  = 


12n  <j>; 


(4.295) 


(4.296) 


As  in  previous  paragraphs,  using  these  large  sample  results,  confidence  in¬ 
tervals  can  be  constructed  for  the  various  parameters. 

In  Tal's  paper  (Reference  32),  least  squares  estimates  of  D  and  <j>  are  derived 
from  the  times  between  error  occurrences  (X^'s)  and  the  times  of  error  occurrences 

(t/s).  Specifically,  for  the  X^'s,  minimize  the  following: 

2 


Si(D,40  =  £  (x  -  -4tt 

i=l\x  D4)1  1 


(4.297) 


Taking  the  partial  derivatives  of  Si  with  respect  to  D  and  <(>,  these  are  obtained: 


9Si  /  n  X.  n  j  N 


(4.298) 


V1  ~ 1}  f  a  -  ,n 
*♦  U=1  D*1  1=1  DH21"1 


(4.299) 
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and 


B.  =  £ 

j=l 


*  LS,G 


(4.306) 


the  laast  squares  estimates  are  the  solutions  of: 
Estimates  -  Least  Squares 

n  n 


LS,G 


S  t.  B.  -  £  C.  B.  =  0 

i=l  11  i=l  1  1 


(4.307) 


and 


n 


n 


D 


LS,G 


£  t.  C.  -  £  C.2  = 

i=l  11  1=1  1 


=  0  . 


(4.308) 


Using  the  results  of  Paragraph  4.2.3,  the  asymptotic  variances  of  these  vari¬ 
ous  estimates  can  be  established.  For  large  n,  the  estimates  have  a  covariance 
matrix  of  the  form: 
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(4.310) 
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4 . 2 . 6 . 1  Modified  Geometric  "De-Eutrophication"  Model .  The  only  extension 
of  the  Geometric  Model  that  is  considered  is  due  to  Lipow38  and  discussed  in 
Sukert.17  The  extension  is  made  to  relax  the  assumption  of  an  infinite  number  of 
errors  teing  presert  in  the  code.  The  model  assumptions  are: 

Model  Assumptions 

(a)  All  errors  do  not  have  the  same  chance  of  detection, 

(b)  The  detections  of  errors  are  independent, 

I 

(c)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage,  and 

(d)  The  error  detection  rate  during  the  ith  time  interval  of  testing  is: 

n'  1 

Z(t)  =  Dtp  1  for  t.^  <  t  <  ti  (4.320) 


where  D  and  0  <  tj>  <  1  are  as  in  the  previous  paragraph  and  n^_1  is  the  cumulative 

number  of  er  ors  found  up  to  the  ith  interval  of  testing.  Ths  form  that  this 
hazard  rate  f  ’notion  takes  is  given  in  Figure  4-4. 
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FIURE  4-4.  MODIFIED  GEOMETRIC  "DE-EUTROPHI CATION"  PROCESS 
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The  development  of  the  MLEs  for  $  and  D  proceeds  in  essentially  the  same 
manner  as  for  the  Geometric  Model.  The  resulting  estimates  are  obtained  as  the 
solutions  of  the  following  pair  of  equations: 


Estimates  -  Maximum  Likelihood 


n 


MG 


m 

^  n  *  m 

£  *£\ 

i=l 


(4.321) 


and 


MG 


m  a  m  a 

£,  “i-1  =  dmg  £  Vi  ♦mg 

1=1  1=1 


vr1 


X. 

i 


(4.322) 


where  m  is  the  number  of  testing  intervals  of  each  length  X.,  and 

m 

n  =  <Pn.  is  the  total  number  of  errors  discovered.  Notice  that  the  hazard  rate 

i=l  1 

function  and  the  estimates  become  those  of  the  previous  section  when  n^_j  =  i-1. 
The  MLEs  of  the  MTBF  and  the  reliability  of  the  program  after  m  intervals  of 
testing  are: 


Reliability  =  R(t)  =  xpJ 


and 


a  /a  n 

"DMG  TMG  ^ 


t  <t 
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MTBF  = 


k* 


A  A  tl_ 

dmg*mg 


(4.323) 

(4.324) 
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The  estimated  degree  of  "purification"  for  a  program  is  obtained  as  in  the  pre¬ 
vious  paragraph  as: 


a  a  n 


PL  = 


Z(to)  -  Z(t„) 
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4.2.7  Geometric  Poisson  Model 


This  model  was  also  proposed  by  Moranda41  as  an  alternative  to  the  Geometric 
Model  if  the  reporting  of  the  software  error  detections  is  on  a  periodic  basis. 
As  in  the  previous  Modified  Geometric  Model,  only  the  numbers  of  error  occurr¬ 
ences  per  testing  interval  are  needed.  Unlike  the  previous  model,  however, 
the  testing  intervals  are  all  assumed  to  be  the  same  length;  e.g.,  a  testing 
period  is  composed  of  a  day,  week,  etc.  Additionally,  since  the  model  assumes  a 
constant  rate  of  error  occurrence  during  a  time  period,  the  model  is  best  applied 
to  situations  in  which  the  length  of  the  reporting  period  is  small  in  relation¬ 
ship  to  the  overall  length  of  the  testing  time.  The  model  assumptions  are: 

Model  Assumptions 


(a)  There  is  a  nonfinite  number  of  errors. 


(b)  The  detections  of  errors  are  independent. 

(c)  The  errors  do  not  have  the  same  chance  of  detection. 

(d)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(e)  During  the  ith  time  period,  the  number  of  errors  detected,  f^,  during 

that  period  follows  a  Poisson  distribution  with  parameter  *  where  D  is  the 
initial  detection  rate  and  ^  is  the  constant  of  proportionality  where  0  <  <)>  <  1. 

(f)  Each  error  discovered  is  either  corrected  or  not  counted  again. 

From  assumption  (e),  the  detection  rate  follows  a  geometric  progression  from 
one  testing  period  to  the  next.  Initially,  the  detection  rate  is  ~  constant  D. 
After  the  first  reporting  period,  the  detection  rate  is  assumed  proportional  to 
the  initial  rate,  i.e.,  it  is  then  $D,  and  so  on.  The  hazard  rate  for  this 
model  is: 


Z(t)  =  DcJ)1 


(4.332) 


for  t^_j  <  t  <  t^  during  the  ith  time  period.  Notice  how  this  compares  to  the 

hazard  rate  functions  for  the  Geometric  and  Modified  Geometric  Models.  Here  the 
t^'s  are  fixed,  while  for  the  Geometric  they  were  random. 

Since  the  number  of  error  detections  in  a  reporting  period  follows  a  Poisson 
distribution,  the  likelihood  function  for  the  m  reporting  periods  is: 


i-lW 


L(f i , • • . »f  )  =  n 

m  i=l 


(4.333) 
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Large  sample  variances  of  these  estimates  can  be  calculated  in  the  usual  way 
and  are  found  to  be: 


and 


r-  i  dgp  T"'1  '  i+i  *  ..1 

{m^-LS1^  ♦  £  « -  »<*  -  »*ij  • 


var 


var  J  <j> 


GPf  ~ 


'GP 
m 

z  ;i-» 

i=l  9GP 

- K - 


AD 


GP 


cov(<PGp>DGp) 


m-1 

-  £  i  ii;1 


i=0 


GP 


(4.341) 


(4.342) 


(4.343) 


l-S 
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«2S3 
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r.i 

fi 


Ig 
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where 


m  ^ ,  t-  .*  a  , 

V  1  m-1  i 

„  £>  -  ?  i  *c 

-  -  i=0 

n  mm 


i+1 

GP 


m 


(b^ 

yGP 


+  E  d  "  DU  “  2)<J>; 

i=l 


»]  ■  it 


i  <() 


i-1 


The  least  squares  estimates  are  obtained  by  minimizing: 
S(<(»,D)  =  E  (f.  -  D^1”1)2  . 


i=l 


(4.344) 


(4.345) 


Taking  partial  derivatives  and  setting  them  equal  to  zero,  the  least  squares  esti¬ 
mates  are  obtained  as  the  solutions  of: 


Estimates  -  Least  Squares 

m  *. .  , 

A  5,  fi*GP,LS 

x=l  ’ 


D 


GP’LS  f'.  2(i-l) 

^  GP,LS 
1=1  ’ 


L  fi^Gl%LS  )(l  -  iop  LS) 

izl - - (A - ZLi±PJ  (4.346) 

{,  22m- 1  \ 

\^~^GP,LS/ 


and 


m  A  ■  ^  u  A  IP  Aq  /  '  a  1 

E  (i  -  l)fi  4>GP  LS  -  DGP  LS  £  C1  “  1  ^ ^GP , LS 
1=1  ’  *  1=1  ’ 


j.(i-l)  =  o, 


(4.347) 


4-65 


I  (.'• 


r./j 

C'.y 

*  *  < 

'  •  *  * 

I* " 


'••AitM/TU*  a.ZL+JlLyijSItLPJiJ'JiM''*  >7  f*  *M»  ^XmSi \*h4 *W*  '%W«^£*i  V»  TMWJW 


1.IU 


^GP , LS^GP , LS 
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(3 


Large  sample  variances  can  be  developed  for  these  estimates  as  done  in  pre¬ 
vious  paragraphs.  The  estimate  of  the  expected  number  of  error  detections  in  the 
(m  +  !)st  time  interval  is  established  as  either: 


Expected  Number  of  Errors  in  Interval  m+l  -  Dgp^Qp 


or 


(4.348) 


1 


4.2.8  Schneidewind1  s  Model 

Norman  Schneidewind43  proposed  a  generalized  model  which  includes  the  Geo¬ 
metric  Poisson  as  a  special  case.  The  basic  philosophy  of  this  model  is  that  as 
the  testing  progresses  over  time,  the  error  detection  process  changes  and  hence, 
recent  error  counts  are  usually  of  more  use  than  earlier  counts  in  predicting 
future  error  counts.  Three  approaches  are  employed  in  utilizing  the  error  count 
data.  Suppose  there  are  m  intervals  of  testing  and  f.  errors  were  detected  in 

the  ith  interval,  one  of  the  following  can  be  done.  1 

(a)  Utilize  all  of  the  error  counts  for  the  m  intervals. 

(b)  Ignore  the  error  counts  completely  from  the  first  s  -  1  time  intervals 
(2  <  s  <  m) ,  and  only  use  the  data  from  intervals  s  through  m. 

(c)  Use  the  cumulative  error  count  from  intervals  1  through  s  -  l,i.e., 
s-1 

W ,  ,  -  £  f.  and  the  individual  errors  counts  from  interval  s  through  m. 

s-i  .=1  i 

Schneidewind  argues  that  approach  number  1  is  applicable  when  one  feels  that  the 
error  counts  from  all  of  the  intervals  are  useful  in  predicting  future  counts. 
Approach  number  2  is  to  be  used  when  it  is  felt  that  a  significant  change  in  the 

error  detection  process  has  occurred  and  thus  only  the  last  m  -  s  +  1  intervals 

are  useful  in  future  error  prediction.  The  last  approach  is  an  intermediate  one 
between  the  two  others.  Here  it  is  felt  that  the  combined  error  count  from  the 
first  s-1  intervals  and  the  individual  counts  from  the  remaining  are  representa¬ 
tive  of  the  error  detection  behavior  for  future  testing  intervals.  The  model 
assumptions  are; 

Model  Assumptions 

(a)  The  number  of  errors  detected  in  one  interval  is  independent  of  the 
error  count  in  another. 


I 


Si 


(b)  The  error  correction  rate  is  proportional  to  the  number  of  errors  to  be 
corrected. 
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(c)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 


(d)  The  mean  number  of  detected  errors  decreases  from  one  interval  to  the 


next. 

(e)  The  intervals  are  all  of  the  same  length. 

(f)  The  rate  of  error  detection  is  proportional  to  the  number  of  errors 
within  the  program  at  the  time  of  test.  The  error  detection  process  is  assumed  to 
be  a  nonhomogeneous  Poisson  process  with  an  exponentially  decreasing  error  detec¬ 
tion  rate.  The  rate  of  change  is  taken  to  be  of  the  form 


dA  =  aexp  {-Pi}  ( 

for  the  ith  interval  where  a  >  o  and  p  >  o  are  the  constants  of  the  model. 
From  assumption  (f),  the  cumulative  mean  number  of  errors  is  therefore 
Di  =  |  t1  "  exp{-pi}]  ,  ( 

so  that  for  the  ith  interval,  the  mean  number  of  errors  is 

mi  =  D.  -  Di_1  =  |  [exp(-p(i  -  1))  -  exp(-pi)].  I 


(4.350) 


(4.351) 


(4.352) 


The  likelihood  function,  assuming  a  Poisson  process,  is  then  developed  as 


M  s,^  exp(-M  m  m.  "'"exp  (-m. 

L(f  f  _ -  i— n  1— -i. 

Mtr***'V  F  . !  f.! 

s-1  i=s  1 


(4.353) 


where  Ms-1  is  the  mean  number  of  errors  in  the  interval  1  through  s-1  with  s 
chosen  as  an  integer  value  in  the  range  2  <  s  <  m. 

Using  the  fact  thats 


"’i  “  f  [exp(-p(i  -  1))  -  exp(-pi)] 

=  p  exp(-p(i  -  1) [ 1  -  exp(-p) ] 


(4.354) 


„  a 


[exp(0)  -  exp(-(s  -  1)P)] 


[1  -  exp(-(s  -  l)p)]  , 


(4.355) 
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Gephart  et.al.18  established  that  the  MLEs  for  a  and  p  are  then  obtained  as: 
Estimates  -  Maximum  Likelihood  (Approach  c.) 


and 


Ps  =  My) 

m 


a  = 
s 


.£  9 


1  -  exp  | -pmj 

where  y  is  the  solution  of  the  polynomial  equation, 


(s  -  1)F  .  F  mF 

S~1  +  s,m  _  m 


y8’1-! 


y  -  1  m  , 
3  y  -  1 


=  A 


which  can  be  simplified  to: 

? 

s,m 


AyS+m  “  (A  +  m)yS+m_1  -  (A  +  sFs.1  -  I^)/*1 


(4.356) 


(4.357) 


(4.358) 


with 


+  (A 

+  F 

+  sF  , 

s  ,m 

s-1 

(A  + 

F  m  -  mFjy3"1 
s  ,m  m 

-  (A 

+  sF  -  +  F 
s-1  s,m 

m 


(4.359) 


4*\i 

Ki 


m-s 

A  =  53  (s  +  i  -  l)f 

i=0 


m 

f  =  y  f. 

s,m  f-»  i 

>  i=s 


s+i 


(4.360) 

(4.361) 


If  s  is  set  to  1  and  F 


s 


(4.358)  gives: 


m 

=  F  =  53  f .  ,  then  the  polynomial  in  equation 
>ra  m  i=i  1 


F  mF 

m  -  m 

y-1  />_! 


=  A. 


(4.362) 
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a  = 
s 


r  mi  \ 

a=s  / 

- c - „ - 

1  -  ext>{-fSm] 


where  y  is  the  solution  of  the  equation 

Aym-S+2  -  (A  +  iSin,)y"'s+1  *  «•  -  ■  +  I)*.,.  -  *)y 


+  (A  +  F  -  (m  -  s  +  1)F  )  =  0 

s,m  s  ,m 


for  y  >  1  where 


(4.365) 


(4.366) 


m-s 

A  =  E  ifs+i* 

i=0  8+1 


(4.367) 


From  the  MLEs,  various  other  parameters  can  be  estimated  as  seen  in  the  following. 


Expected  Number  of  Errors  in  the  (m  +  l)st  interval  of  testing 


=  mi+l  =  |  [exp{-0i)  "  exp{-p(i  +  1)}] 


(4.368) 


Time  to  detect  a  total  number  of  M  errors 


=  log  [a/ (a  -  pM)]/p  ; 


(4.369) 


This  simplifies  to: 


Aym  -  (A  +  Fm)ym  +  (mFm  -  A)y  +  (A  +  Fffl  -  mFm)  =  0, 


(4.363) 


for  y  >  1. 


The  MLEs  derived  under  these  conditions  are  for  approach  (a)  where  all  error  count 
data  are  used. 

In  equation  (4.363),  if  m  -  s  +  1  is  substituted  for  m  and  the  subscript  of 
the  f  Js  is  modified  in  the  expression  for  the  summations  to  make  f  the  first 
error  count,  the  MLEs  for  approach  (b)  are  obtained  where  the  first  s  -  1  error 
counts  are  ignored.  For  this  case,  the  MLEs  are: 


Ps  =  to(y) 


(4.364) 
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the  correction  rate  for  the  ith  interval 
=  a  exp  {-p(i  -  Ai)}  ; 


(4.370) 


where  Ai  is  the  lag  time  between  the  detection  of  errors  and  their  correction, 

i.e,  the  time  to  correct  D.  -  C.  errors  where  C.  is  the  cumulative  number  of 
’  1  ’  ■> 


C.  errors  where  C.  is  the  cumulative  number  of 
i  i 


errors  corrected  up  through  the  ith  interval. 

All  of  these  parameters  are  estimated  by  substituting  the  appropriate  MLEs 
for  a  and  p.  If  the  lag  Ai  is  unknown,  it  can  be  estimated  by  finding  a  value 
for  Ai  such  that 


C.  =  D.  .. 
x  i-Ai 


,  i  >  Ai 


(4.371) 


using  the  empirical  data. 


If  approach  (b)  or  (c)  is  used,  a  determination  for  s  needs  to  be  made. 
Schneidewind  suggests  letting  s  =  2,...,m  and  finding  the  MLEs  for  each  value  of 
s.  For  each  pair  of  estimates,  the  computed  sum  of  weighted  squared  deviations 
between  the  error  estimates  m^  and  the  observed  counts  f^  for  all  i  is  computed 

and  the  one  yielding  the  smallest  sum  is  the  chosen  s.  The  weighted  sum  is  given 


SDW  =  ]C  exp  (pi)  j|  exp(-pi)J  jexp(p)  -  1  j  -f  ±  J 


(4.372) 


Schneidewind  also  suggests  that  to  decide  among  which  of  the  three  approaches  to 
use  [(a),  (b),  or  (c)],  the  unweighted  sum  of  squares 


=  £  [i 

i=m+l  LP 


exp(-pi)  {exp(P)  -  1}  -f. 


(4.373) 


is  computed  for  each  approach  (i.e.,  a  and  p  are  replaced  by  their  estimates  for 
the  respective  model).  M  is  some  specified  future  time.  The  unweighted  sum  of 
squares  is  calculated  between  the  observed  counts  and  the  expected  counts  over  the 
next  M  -  m  intervals.  The  approach  yielding  the  smallest  sum  and  hence,  yielding 
the  smallest  differences  between  predicted  and  actual  values  is  the  one  chosen. 

Gephart  et.al.18  show  that  the  models  under  approaches  (a)  and  (c)  (with 
s  =  2)  are  equivalent  to  the  Geometric  Poisson  of  Paragraph  4.2.7.  If 


(4.374) 
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P  -  -£n<J>  , 


(4.375) 


and  these  are  substituted  into  Schneidewind' s  Model,  it  becomes  the  Geometric 
Poisson  Model  where  D  and  <j>  are  the  parameters  of  that  model.  In  an  equivalent 
manner,  if 


B  »  |  [1  -  e'1*] 


(4.376) 


(j)  =  e 


(4.377) 


and  they  are  substituted  into  the  Geometric  Poisson,  it  becomes  Schneidewind' s 
Model. 

The  data  required  to  implement  any  one  of  three  models  are: 

Data*  Requirement 

The  error  counts  for  each  of  the  m  intervals  of  testing. 


4.2.9  Nonhomogeneous  Poisson  Process 

The  Nonhomogeneous  Poisson  Process  (NHPP)  Model  was  proposed  by  Amrit  Goel 
and  Kazu  Okumoto  (References  44,  45,  and  46).  Following  other  models  that  have 
been  considered  (see  Paragraphs  4. 2. 3. 2,  4.2.5,  4.2.7,  and  4.2.8),  this  model  as¬ 
sumes  that  the  error  counts  over  nonoverlapping  time  intervals  follow  a  Poisson 
distribution.  The  expected  number  of  errors  for  the  Poisson  process  in  an 
interval  of  time  is  assumed  proportional  to  the  remaining  number  of  errors  in  the 
program  at  that  time.  Specifically,  the  model  assumptions  are  as  seen  in  the 
following. 


Model  Assumptions 

(a)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(b)  The  numbers  of  errors,  (f i,f2» • *  * *fm) >  detected  in  each  of  the  respec¬ 
tive  time  intervals  [ (0,ti) ,  (tx ,t2) , (t2 ,t3) , . . . , (ti_1,t.) , . .  •  .(t^ptJ  J  are  inde¬ 
pendent  for  any  finite  collection  of  times  ti  <  t2  <,..., <t  . 

(c)  Every  error  has  the  same  chance  of  being  detected  and  is  of  the  same 
severity  as  any  other  error. 


’>  v  v  v*  v  v  v-  v  v v  v  v  v  •  >  v  ■“.*  ^  • .  ■ .  v.v.  v ■ .  * .. •. 
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(d)  The  cumulative  number  of  errors  detected  at  any  time  t,  (N(t)),  fol¬ 
lows  a  Poisson  distribution  with  mean  m(t) .  The  mean  m(t)  is  such  that  the  ex¬ 
pected  number  of  error  occurrences  for  any  time  (t,t  +  At)  is  proportional  to  the 
expected  number  of  undetected  errors  at  time  t. 

(e)  The  expected  cumulative  number  of  errors  function,  m(t),  is  assumed  to 
be  a  bounded,  nondecreasing  function  of  t  with 

m(t)  =0  t  =  0 

m(t)  =  a  t-x* 

where  a  is  the  expected  total  number  of  errors  to  be  eventually  detected  in  the 
testing  process. 

Note  that  f.  =  N(t.)  -  N(t.  ,).  The  NHPP  differs  from  some  of  the  other 
1  v  i  l-l 

Poisson  Models  considered  in  that  this  model  treats  the  initial  error  content  of 
a  program  as  a  random  variable  while  some  of  the  others  assume  it  is  a  fixed 
constant.  Also  the  time  between  the  (i  -  l)st  failure  and  the  ith  failure  de¬ 
pends  upon  the  time  to  failure  of  the  (i  -  l)st  rather  than  being  independent 
of  it. 


From  assumptions  (d)  and  (e),  for  any  time  period  (t,t  +  At) 

m(t  +  At)  -  m(t)  »  b{a  -  m(t)}At  +  O(At)  (4.37T.) 

where  b  is  the  constant  of  proportionality  and  0  as  At  0.  By  letting 

At  ->  0,  the  mean  function  satisfies  the  differential  equation 

m'(t)  =  ab  -  bm(t) .  (4.379) 

Under  the  initial  condition  m(0)  =  0,  the  mean  function  is 

m(t)  =  a(l  -  e“bt)  .  (4.380) 


Thus , 

Pr{N(t)  =  n} 


n! 


with 

m(t)  =  a(l  -  e"bt)  .  (4.381) 

For  f ^  =  N(t^)  -  N(t^)  and  the  error  counts  being  independent,  the  likelihood 
function  is  therefore: 


L(f 


,f  )  = 

’  in 


m 

n 

i=l 


[>*n(l-i) 


m 


(ti_1)]  1exp{ra(ti-1)  -  m(ti)} 


f .! 

l 


(4.382) 
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■bt .  f . 


-bt.  -bt. 


■e  1)]  'SxpJaCe  1-  e 


m  m 

AnL(flt . . .  ,fm)  =  £  f^na  +  £  f.£n(e 
i=l  i=l 


-bt.  .  -bt. 
"-1  -e  X) 


3  :• 

m  -bt .  -bt .  , 

m 

p  '■ 

+  a  Z  («  1  -e  1_1)  “ 

E 

R 
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i=l 

i=l 
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R  , 
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m 

i fe 

8£nL  _  6i  1  ,,  ~btms 

9a  “  a  -  0  “  e  ) 

(4.383) 
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Thus,  the  MLEs  are  the  solutions  to  the  following  system  of  equations: 

Estimates  -  Maximum  Likelihood 

m 


£  £i 

i=l 


NHPP  (1  -  e  "bNHPPtm) 


r  .  /  -b  t.  -b  t.  J 

,  "DNHPPm  ,  L  NHPP  1  -t.  te  NHPP  1 

t  e  }.  f .  f .  I  t.e  l-l  , 

m  j=i  1  _  f.  1  \  1 _ l 

(l-e  "bNHPPtm)  1=1  e"bNHPPti-l  -  e"bNHPPti 
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and 


“bNHPP,LSti-l 


t  h  \tt.  -  tw. 

>ls  £  [(f'bNHPP,LSti-l  .  /bNHPP,LSti-t^tie-bNHPP,LSti  .  ^  ^NHPP.Lsh-jj 


In  the  report  by  Schafer38,  the  large  sample  variances  and  covariances  are 
derived  for  the  MLEs  and  the  least  squares  estimates.*  Provided  a  is  large,  the 
variances  and  the  covariances  of  the  MLEs  are: 


var{a} 


MiPP 


var{b} 


A  ** 

(t.  "  t_1)2exp(-b(ti  +  t._1))  -bt 


bt.  .  -bt 
e  -e  1 


0 


t2e  m 
m 


and 


NHPP 


A  ~  1  "bt 

cov(aNHPP,bNHPP)  ~  5  (-tme  *  » 


(4.396) 


(4.397) 


(4.398) 


with 


-bt  /  m 
A  =  (1  -  e  m)  [  Z 

i=l 


-b(t.  +  t.  .) 

(t  t.  ,)2e  1  1-1  -bt  \  -2bt 

1  1-1  .  t2p  ml  .z  m 

me  /  V  (4.399) 


(■ 


-bt .  .  -bt . 

i-l  i 

e  -  e 


For  the  least  squares  estimates,  the  large  sample  variances  and  covariances  are 
given  in  the  covariance  matrix: 


<:  ~i_  £ 

°NHPP,LS  “  A2  1 


*  The  author  would  like  to  express  his  thanks  to  Mr.  Vijaya  K.  Srivastava  who 
pointed  out  the  errors  in  the  original  formulas  and  provided  the  correct  ones. 
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where 


*NHPP,LS 

"’NHPPjLS  "\r  1  ’ 

'DNHPP,LS, 


m 


A  =  £  a  1  e 

i=l 


-bti\3 


-  e 


m 


B  =  X]  a‘l  e 


i-1 


-bt 


C  =  £  a-  V  e 


i-1 


i=l 
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(4.402) 

(4.403) 
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(4.405) 


and 


.  r  /  -bti  -bti-t  -  /  -bti  -btt-v-bti-i  -b* 

E  [•(v  -  t..,e  j|  -  E  A4!6  ■  h-1'  /\e 


with: 


.  I  -bti  'bti-l)  (  -bti-l  -bti)  A  (  'bti-l  'bti)  2 

E  a^.e  -t,,.  )\e  -  '  )  £  (.  -•  J 

1  (4.406) 

This  large  sample  result  is  derived  utilizing  the  results  of  Paragraph  4.2.3 


-bti-l  -bti! 


gj  =  «  \e  -  e 


and 


’bti-l  "bti\ 


o?  =  a  |e  -  e 


(4.407) 

(4.408) 
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Equations  (4.396  through  4.406)  can  then  be  used  to  construct  large  sample 
confidence  intervals  for  the  parameters  by  replacing  any  unknown  in  the  variances 
with  their  respective  estimates. 


Goel  and  Okumoto44  derive  MLEs  for  a  and  b  based  upon  the  individual  times  of 
error  occurrences.  If  represents  the  failure  time  of  the  ith  error,  they  show 

in  their  report  that  the  MLEs  are  the  solutions  to  the  following  system  of  equa¬ 
tions: 


A 


n 

- 


1  -  exp(-b(J0sn) 


(4.409) 


and 


=  £  sk  +  aGOV*p  ('bGOsn>  (4-A10) 

dgo 

where  n  is  the  total  number  of  errors  detected. 


Using  this  formulation,  Goel  and  Okumoto  establish  that  if  S  =  s  is  the  time 

n 

of  the  last  failure,  then  the  conditional  reliability  function  of  X  (the  time 
between  the  nth  and  (n  +  l)st  failures)  is  given  by: 

Xn  U  |S"  ’  S>  =  *  lX“+1  J  *  lS»  =  81  <4'4n) 

=  P  {software  is  operational  for  at  least  x  amount 
of  time  given  s  amount  of  testing}  (4.412) 

=  exp  [-a{e  bs  -  e  k(s+50}]  (4.413) 

Okumoto  and  Goel45  utilize  this  reliability  to  determine  an  optimal  release 
time  for  a  software  program.  Testing  can  continue  until  the  desired  reliability 
R  =  Ry  (x  IS  =  s)  is  achieved  for  a  specified  operational  time  of  x  or  the 
n+1  n 

required  testing  time  s  can  be  determined  for  a  desired  reliability  for  a  speci¬ 
fied  operational  time.  If 

R  =  exp[-a{e”bs  -  e"b(s+x)}]  (4.4l4) 

=  exp[-m(x)e  bs]  ,  (4.415) 


with 


m(x)  =  a(l  -  e  bx) 


then  the  desired  testing  time  to  achieve  the  specified  R  and  x  is: 


8  =  b 


£n^rn(x)^  -  ^£n  ^£n  Qj 


(4.416) 

(4.417) 
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Estimates  based  upon  previous  testing  data  are  used  for  a  and  b. 

In  their  paper,  Okumoto  and  Goel  also  determine  an  optimum  release  time  based 
upon  cost  considerations.  Suppose: 

Ci  =  cost  of  fixing  an  error  during  testing, 

C2  =  cost  of  fixing  an  err.'r  during  operational  use  (C2  >  Ci)> 

C3  =  cost  of  testing  per  unit  time, 

t  =  software  life  cycle  length,  and 

T  =  software  release  time  for  testing. 

Since  m(t)  is  the  cumulative  expected  number  of  errors  in  the  interval  (0,t),  then 
the  total  expected  cost  is: 

C(T)  a  Cim(T)  +  C2  [m(t)  -  m(S)]  +  C3T.  (4.418) 

Differentiating  the  expression  with  respect  to  T,  then 

C'(T)  =  Cim'Cn  -  C2m'(T)  +  C3  (4.419) 


where 


ra'(T)  =  abe“bT  . 


(4.419) 


(4.420) 


Setting  the  right-hand  side  of  this  equation  equal  to  zero,  the  following  is 
obtained: 


-bT  _ 
&  — 


C2  -  Ci 


(4.421) 


Okumoto  and  Goel  establish  in  their  paper  (Reference  45)  that: 

C3 

(a)  If  ab  >  ;; - — ,  then  there  exists  a  unique  feasible  solution  to  equa- 

C2  -  Ci 

tion  (4.421)  and  the  optimum  release  time  is: 

T*  =  min{T0,t}  (4.422. 


where 


To  =  i  to 


ab(C2-Ci) 


(4.423 


while 


If  ab  < 


(4.424 
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then 

T*  =  0.  (4.425) 

The  last  comment  concerning  this  model  is  that  if  the  testing  intervals  are 

all  o£  the  same  length,  say  T,  then  this  model  is  equivalent  to  the  Geometric- 

Poisson  Model  and  hence,  Schneidewind' s  Model  with  s  =  1  of  Paragraphs  4.2.7  and 

4.2.8  respectively.  If  all  of  the  testing  intervals  are  of  the  same  length,  then 

the  time  of  the  ending  of  the  ith  testing  interval  is  t.  =  iT.  The  joint  density 

of  the  f . ' s  then  becomes :  1 

i 


£(f> . V 


£ 

m  [m(ti)  -  m(ti_1)"|  1exp  jmCt^)  -  m(ti)} 

n  r  I  — ^ 


i=l 


f.  ! 
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.  [ 
=  n  L 
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/  l-l  i 

a  (e _ -  e 


l)  1exp  {  a(< 
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i  l 

e  -  e 


i=l 


f .  f 
i 


(4.426) 


(4.427) 
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(4.428) 


m  I"  -bTv  -b(i-l)Tl  [  -bTv  -b(i-l)T 

-  n  la(1"e  )e  J  exp  lra(l  -  e  )e 
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(4.429) 


Notice  that  if: 


and 


b  = 

T 


(4.430) 


(4.431) 


where  D  and  $  are  defined  in  Paragraph  4.2.7,  then  the  joint  density  function 
becomes : 


m  n,ii“  1  ^i  I 

n  [frfr  exp  I-D4>  J 

11  r  i 


i-1 


f . 

i 


(4.432) 


the  Geometric  Poisson.  Likewise,  since  the  Geometric  Poisson  and  Schneidewind' s 
Model  with  s  =  1  are  equivalent,  utilizing  the  relationships  established  between 
the  two  models  in  the  last  paragraph, 


and 


a 

3  =  P 


b  =  T 


(4.433) 

(4.434) 
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are  found  to  be  the  relationship  between  Schneidewind ' s  and  the  NHPP  Models. 

The  data  required  to  implement  this  are: 

Data  Requirements 

(a)  The  error  counts  in  each  of  the  testing  intervals,  (i.e.,  the  f/s). 

(b)  The  times  the  testing  intervals  end,  (i.e.,  the  t^'s). 

(c)  The  time  of  error  occurrences  (i.e.,  the  s.'s)  if  an  optimal  release 

time  is  desired.  1 


4.2.10  Duane's  Model 

The  next  model  considered  also  employs  a  nonhomogeneous  Poisson  process  for 
the  error  counts.  This  model  was  originally  proposed  by  J.  T.  Duane4'  as  a  hard¬ 
ware  reliability  growth  model.  Duane  observed  that  the  cumulative  failure  rate 
versus  cumulative  testing  time  when  plotted  on  An-An  paper  tended  to  follow  a 
straight  line  for  a  number  of  systems  developed  at  General  Electric.  This  model 
has  been  applied  with  some  succebs  to  software  reliability  modeling  by  Evaluation 
Associates,  Inc.,  (References  48  through  50).  The  specific  assumptions  for  this 
model  are  given  in  the  following. 


Model  Assumptions 

(a)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage, 

(b)  Every  error  has  the  same  chance  of  being  detected  and  is  of  the  same 
severity  as  any  other  error, 

(c)  The  error  occurrences  are  independent,  and 

(d)  The  cumulative  number  of  errors  detected  at  any  time  t,  [N(t) ] ,  follows 
a  Poisson  distribution  with  mean  m(t).  The  mean  function  is  taken  to  be  of  the 

form  m(t)  =  \t^. 

From  the  assumptions,  it  can  be  seen  that  if 

m(t)  _  At^  _  Expected  number  of  errors  by  time  t  n 

~t  ~  t  total  testing  time  (4.435) 


is  plotted  on  An  -  An  paper  versus  time,  or  conversely,  if 


Y  = 


=  An 


=  AnA.  +  (8  -  l)An  t 


(4.436) 
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is  plotted  on  regular  paper  versus  2n(t),  a  linear  relationship  relating  the  two 
is  obtained.  That  is, 

Y  =  a  +  bX  (4.437) 

with  a  =  JinA,  b  =  P  -  1,  and  X  =  £n(t)  for  the  latter  case. 

The  rate  at  which  errors  are  occurring  is: 

(4.438) 


dm(t)  _  XfttP-l 
dt  -  Apt 


Hence,  the  MTBF  is 


1 


Apt 


P-1 


Also 


if  p  >  1,  there  is  no  improvement  in  the  software  as  time  progresses.  Crow51 
shows  the  MLEs  for  A  and  p  are: 


Estimates  -  Maximum  Likelihood 


= 


_  n 


(4.439) 


and 


PD  ~  n-] 


n 


(4.440) 


E  2n(t  /t.) 
i=l  n  1 

where  the  t^'s  are  the  observed  failure  times  and  n  is  the  number  of  software 
errors  detected. 

The  MLE  of  the  MTBF  for  the  (m  +  l)st  error  occurrence  is  then: 

HW  =[^PDtnPD'1  ]  'l  =  V  • 

nPD 

Crow51  also  provides  a  table  of 
P  |  MTBF/MTBF  <  CQ  |  =  a 


(4.441) 


(4  442) 


which  can  be  used  to  construct  a  100X(1  -  a)  percent  confidence  interval  for  the 
MTBF. 
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Least  squares  estimates  for  a  =  JdnA  and  b  =  p  -  1  for  the  equation 

/,_  (Expected  number  of  errors  by  time  t  \ 

An  V  - total  testing  time  - /  "  a  +  bAn(t) 


can  be  achieved  in  the  standard  manner  as: 


(A. 443) 


Estimates  -  Least  Squares 


aD,LS  =  Y  "  bD,LSX 


(4.444) 


D,LS 


where 


n  n  n 

Z  x  Y  -  £  X  £  Y 

i=l  1  1  i=l  1  i=l  1 
•  £  1  -  (£'.)• 


(4.445) 


Xi  - 


^  ■  >ik) ' 


(4.446) 


(4.447) 


Various  confidence  intervals  for  the  parameters  of  the  linear  model  can  be 
constructed  in  the  usual  way. 

The  only  data  required  to  implement  this  model  are: 


Data  Requirement 


The  times  of  error  occurrences. 


4.2.11  Execution  Time  Model 

The  next  model  considred  is  one  that  has  been  applied  to  the  greatest  number 
of  software  development  programs.  This  is  a  model  developed  by  John  Musa  of 

Rail  T.aWgtnrias  52  *  53 » 54 1  55 » 56  Tho  ^nf»r«st'ina  a»n»rt  ftf  fhis  mnH»l  is  that-  if 


Bell  Laboratories.^*63*5** 55,56  The  interesting  aspect  of  this  model  is  that  it 
is  based  upon  the  amount  of  CPU  time  involved  in  testing  rather  than  on  calendar 
(wall  clock)  time;  but,  the  model  attempts  to  relate  the  two.  By  doing  this, 
Musa  is  able  to  model  the  amount  of  limiting  resources  (failure  identification 
personnel,  failure  correction  personnel,  and  computer  time)  that  may  come  into 
play  during  various  time  segments  of  testing.  In  addition,  this  model  eliminates 
the  need  for  developing  an  ei  correction  model  since  the  error  correction  rate 
is  directly  related  to  the  ins*._  taneous  failure  rate  during  testing.  The  spe¬ 
cific  assumptions  for  this  model  are  given  in  the  following. 
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Model  Assumptions 


(a)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(b)  The  detections  of  errors  are  independent. 

(c)  All  software  failures  are  observed. 

(d)  The  execution  times  between  failures  are  piece-wise  exponentially 
distributed  (i.e.,  the  hazard  rate  is  a  constant  that  changes  only  at  each  error 
correction) . 

(e)  The  hazard  rate  is  proportional  to  the  number  of  errors  remaining  in 
the  program. 

(f)  The  fault  correction  rate  is  proportional  to  the  failure  occurrence 

rate. 

(g)  The  quantities  of  the  resources  (failure-identification  personnel,  fail¬ 
ure  correction  personnel,  and  computer  times)  that  are  available  are  constant  over 
a  testing  segment. 

(h)  Resource  expenditures  for  the  kth  resource,  Axk>  associated  with  a 
change  in  MTBF  from  T*  to  T2  can  be  approximated  by: 

Axk  ~  0R  At  +  Mk  Am  (4.448) 

where  At  is  the  increment  of  execution  time,  Am  is  the  increment  of  failures  ex¬ 
perienced,  0^  is  an  execution  time  coefficient  of  resource  expenditure,  and 

is  a  failure  coefficient  of  resource  expenditure. 

(i)  Failure-identification  personnel  can  be  fully  utilized  and  computer 
utilization  is  constant. 

(j)  Failure-correction  personnel  utilization  is  established  by  limitation 
of  error  queue  length  for  any  debugger.  Error  queue  length  is  determined  by 
assuming  that  error  correction  is  a  Poisson  process  and  that  servers  are  randomly 
assigned  in  time. 

Assumptions  (g)  through  (j)  are  needed  if  there  is  interest  in  modeling 
resource  allocation  for  the  testing  segments.  Only  (a)  through  (f)  are  needed 
for  reliability  modeling.  In  fact,  (a)  through  (e)  are  assumptions  which  are 
incorporated  into  many  of  the  models  presented  in  this  report.  Later  in  Para¬ 
graph  4.2.11,  an  equivalence  relationship  between  this  model  and  the  Jelinski- 
Moranda  Model  of  Paragraph  4.2.3  is  established. 

Suppose  there  is  an  initial  number  of  N  errors  present  in  the  program. 
Suppose  n  errors  have  been  corrected  after  t  amount  of  testing  (based  upon  CPU 


time)  has  elapsed.  Then  £rom  assumption  (e) ,  the  hazard  rate  function  at  time  t 
is  of  the  form: 

Z(x)  =  fK(N  -  n)  ,  (4.449) 

where  f  is  taken  as  the  linear  execution  frequency  (average  instruction  rate 
divided  by  the  number  of  instructions  in  the  program)  and  K  is  an  error  exposure 
ratio  which  relates  error  exposure  frequency  to  linear  execution  frequency.  The 
error  exposure  ratio  attempts  to  account  for  the  fact  that  code  is  not  executed 
in  a  sequential  manner,  due  to  numerous  loops  and  branches,  and  for  the  variation 
of  the  machine  state.  The  variation  of  the  machine  state  may  cause  an  error 
associated  with  a  particular  instruction  to  be  undetected  on  a  given  execution  of 
the  instruction. 

From  assumption  (f), 

=  BZ(l)  ,  (4.450) 

where  B  is  the  proportionality  constant.  B  is  called  the  error  reduction  factor. 
It  is  the  average  ratio  of  the  rate  of  reduction  of  errors  to  the  rate  of  failure 
occurrence.  Usually  B  is  positive  and  less  than  1  although  there  is  the  situation 
in  which  the  finding  of  the  error  that  led  to  the  failure  of  the  program  leads  to 
the  discovery  of  additional  errors  as  well.  This  creates  a  B  larger  than  1. 

Musa  generalizes  this  relationship  by  considering, 

“V  =  BCZ(t)  (4.451) 


where  B  is  as  before  and  C  is  a  constant  called  the  testing  compression  factor. 
It  is  the  average  ratio  of  rate  of  detection  of  errors  during  testing  to  that 
during  use.  It  attempts  to  account  for  the  greater  stress  that  is  placed  on  a 
program  to  uncover  program  errors  during  the  testing  phase  in  contrast  to  the 
operational  phase.  Usually  C  is  larger  than  1  because  of  this  fact. 

Now  suppose  m  represents  the  number  of  failures  experienced  in  the  process 
of  correcting  n  errors  and  suppose  M  is  the  required  number  of  failures  that  one 
needs  to  experience  to  uncover  all  N  errors  within  the  program.  Then 


n  =  Bm 


and 


N  =  BM 


(4.452) 


(4.453) 


The  previous  equations  can  be  combined  to  obtain: 


dn 

dt 


=  BCZ(t) 


(4.454) 

(4.455) 


=  BC[fK(N  -n)] 


NSWC  TR  82-171 


=  BCfKN  -  BCfKn  , 

i.e. , 

+  BCfKn  =  BCfKN 

or  in  terms  of  the  m's, 

B  —  +  B2CfKm  =  B2CfKM 


or 


dm 

dt 


BCfKm  =  BCfKM  . 


Since  n=m=0att=0,  equation  (4.457)  has  the  solution 
n  =  N  [1  -  exp(-BCfKt) ]  . 
and  equation  (4.459)  has  the  solution 
m  =  M  [1  -  exp(-BCfKl) ]  . 

Since  the  MTBF  is  given  by: 

““■ihr  • 


(4.456) 

(4.457) 

(4.458) 

(4.459) 

(4.460) 

(4.461) 

(4.462) 


it  can  be  reexpressed  as: 

MTBF  =  zk 

1 

■  fK(N  -  n) 

=  fK(K  -  H  +  n‘exp  <-BCfKT»  >  using  e«uati0»  (‘>■‘<<•0), 
1 

"  fKN  exp  (-BCfKt)  ‘ 

If  Tq  is  the  initial  MTBF  when  testing  just  begins,  i.e.,  t  =  0,  then 
T0  =  Initial  MTBF  = 

-  _JL 

~  fKN  ’ 


(4.463) 

(4.464) 

(4.465) 

(4.466) 

(4.467) 

(4.468) 
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Thus , 


MTBF  =  T0  exp  (BCx/NT0)  (4.469) 

for  any  testing  time  T.  As  T  4  »,  the  MTBF  ■+  «  indicating  the  improvement  in  the 
software  as  testing  proceeds. 


The  reliability  of  the  program  at  any  future  time  Xi  given  testing  of  length 
X  is  found  from: 


R(Xi)  =  exp 


[-  Z(x)dx] 
0 


exp  1-XiZCx)) 


(4.470) 


Again  it  can  be  seen  as  x  4  »  causing  both  Xi  and  MTBF  ■>  »,  R(Xi)  4  1  is  obtained. 


From  this  basic  model,  Musa  establishes  some  other  useful  results  (Refer¬ 
ence  52) .  The  number  of  failures  Am  that  must  be  detected  and  corrected  to 
achieve  an  increase  in  MTBF  from  T^  to  Tj  can  be  shown  to  be: 


Am  =  MT0 


(4.471) 


The  additional  execution  time  required  to  achieve  the  increase  is: 


Ax  =  ^  An 


(4.472) 


For  the  implementation  of  this  model  for  a  reliability  analysis,  an  idea 
of  what  the  values  of  these  various  parameters  are  is  needed.  Musa  suggests  that 
initial  estimates  can  be  obtained  from  other  projects  of  a  similar  nature.  For 
some  of  the  parameters  reestimation  can  then  be  made  as  the  testing  progresses. 
The  error  reduction  factor,  B,  can  be  determined  by  taking  data  on  the  number  of 
errors  generated  while  fixing  other  errors.  This  information  could  be  obtained 
from  the  development  of  similar  programs.  Musa  reports  that  B  is  relatively 
stable  for  the  programs  he  considers;  it  is  in  the  range  .94  to  1.00. 

The  testing  compression  factor  C  must  also  be  obtained  in  a  similar  manner. 
If  there  is  no  basis  for  estimation  of  C,  a  conservative  approach  of  taking  C  =  1 
is  advised.  An  initial  value  of  M  can  be  obtained  from  the  relationship  M  =  N/B 
with  N  being  estimated  from  an  idea  of  the  average  error  rate  for  programs  of  a 
similar  nature.  Musa  notes,52  from  a  number  of  other  studies  being  observed,  er¬ 
ror  rates  in  the  range  of  3.36  to  7.98  errors  per  thousand  lines  of  instruction 
with  a  weighted  mean  of  5.43  errors  per  thousand  lines  of  instruction.  In  a  later 
report  (Reference  55),  Musa  employs  an  estimate  of  6.25  errors  per  thousand  lines 
of  instruction.  The  accuracy  of  the  initial  estiamtes  for  N  and  hence  M  do  not 
have  to  be  very  high  since  as  the  testing  progresses,  they  are  reestimated.  The 
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parameter  K,  the  error  exposure  ratio,  must  also  be  estimated  initially  from  pro¬ 
grams  of  a  similar  nature;  however,  like  M  and  N,  it  can  be  reestimated  as  the 
testing  progresses.  Therefore,  the  initial  accuracy  for  this  estimate  can  be  low 

"6 

as  well.  Musa  observes  an  average  value  of  1.31  x  10  for  K  for  the  various  vali' 
dating  projects  he  considered. 

For  the  reestiraation  of  K  and  M,  suppose  Xi,...,X  are  the  times  between 

error  occurrences.  Using  assumption  (d) ,  Musa52  establishes  the  MLEs  for  M  and 
To,  the  initial  MTBF,  as  the  solutions  to  the  following  system  of  equations: 


Estimates 


and 


where 


and 


T0  =  C  X 


M 

ID 


[l  '  l  *] 


M  l 

?  J 

j=M-m+l 


=  <t> 


a  1  m 

*  =  -—  £  (i  -  l»i 

m2  X  i=l  1 


X  = 


m 


ra 


(4.473) 


(4.474) 


(4.475) 


(4.476) 


The  estimate  of  K  is  then  obtained  from  the  relationship: 

t  -  J_  -  JL_  . 

10  “  fKN  "  fKBM  ’ 

i.e. , 


IV  —  - 7V - . 

fBT0M 

If  X  amount  of  testing  has  been  completed  (measured  in  CPU) ,  then  from  the  results 
that  have  been  established  earlier: 


(4.477) 

(4.478) 


MTBF  =  T0  exp  (Ct/MT0)  (4.479) 
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and  Z  „  is  taken  from  a  normal  table  such  that: 
1-- 

2  /  \ 


(4.488) 


Musa's  Model  can  be  related  to  Jelinski-Moranda' s  Model  (therefore,  to  all 
models  which  have  been  shown  equivalent  to  it)  by  letting,  * 


n  =  m 
N  =  M  , 


(4.489) 

(4.490) 


This  is  precisely  the  hazard  rate  for  the  Jelinski-Moranda  Model. 
Also  note  that: 

T0  =  Initial  MTBF  =  ™ 


(4.495) 


which  together  gives  B  =  1, 

(4.491) 

and  fK  =  <p  . 

(4.492) 

If  a  point  X  is 
error,  the  hazard 

chosen  between  the  occurrence  of  the  (i 
rate  function  for  Musa's  Model  is 

-  l)st  error  and  the 

r? 

Z(T)  =  fK(N 

-  (i  -  D) 

(4.493) 

vV*.' 

V-W 

=  4>  (N 

-  (i  -  l))  . 

(4.494) 

m 

becomes : 


(4.496) 


the  correct  expression  under  the  Jelinski-Moranda  Model.  The  MTBF,  after  the 
discovery  of  (i  -  1)  errors  for  Musa's  Model,  was  shown  to  be: 


MTBF  = 


Zlx7  =  fK(N  -  i  +  1) 


t._i  <  T  <  t  . 

With  the  previous  relationships ,  then: 


MTBF  = 


fK(N  -  i  +  I! 


(4.497) 

(4.498) 
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-  1 
~  <j>(N  -  i  +“l) 


(4.499) 


again  the  correct  expression  for  the  MTBF  under  the  Jelinski-Moranda  Model. 


The  importance  of  Musa's  Model  is  in  its  development  of  resource  allocation 
and  the  relationship  between  CPU  time  and  wall  clock  time.  The  resources  (fail¬ 
ure  identification  personnel,  failure  correction  personnel,  and  computer  time) 
influence  the  failure  detection  rate  during  the  testing  process.  At  any  point 
in  the  testing  cycle,  one  of  these  resources  limits  the  other  two  and  thus,  the 
error  detection  rate.  For  example,  if  the  number  of  failure  correction  personnel 
is  insufficient  to  handle  the  errors  detected  by  the  failure  identification  per¬ 
sonnel,  a  backlog  of  errors  develops,  slowing  down  the  testing  process.  Usually 
the  testing  process  involves  from  one  to  three  periods,  each  one  characterized  by 
a  different  limiting  resource.  At  the  start  of  testing,  when  numerous  errors  are 
discovered,  the  limiting  factor  is  the  failure  correction  personnel.  As  the  test¬ 
ing  progresses  and  longer  intervals  between  failures  are  observed,  the  failure 
correction  personnel  utilization  drops,  while  the  failure  identification  personnel 
becomes  the  limiting  factor.  Finally,  at  longer  failure  intervals,  the  use  of  the 
computer  becomes  the  prime  limiting  factor.  Musa's  Model  attempts  to  utilize  the 
knowledge  of  these  limiting  resources  to  relate  execution  time  with  the  passage  of 


calendar  time. 


dtj  dt„  dt^ 

Suppose  2^-,  and  are  the  instantaneous  calendar  time  to 


execution  time  ratios  that  result  from  the  effects  of  each  of  the  resource  con¬ 
straints  taken  alone.  The  index  I  denotes  failure  identification  personnel,  F 
denotes  failure  correction  personnel,  and  C  denotes  computer  use.  An  increment, 
in  calendar  time,  At,  is  taken  to  be  proportional  to  the  average  amount  by  which 
the  limiting  resource  constraints  testing  over  a  given  execution  time  segment; 
that  is, 


At  = 


max 


(4.500) 


From  assumption  (h) ,  the  resource  requirements  associated  with  a  change  in 
MTBF  from  T*  to  T2  can  be  approximated  by: 


Axk  ~  ek  Ax  +  Mk  (4.501) 

where  At  is  the  increment  of  execution  time,  Am  is  the  increment  of  failures  ex¬ 
perienced,  8^  is  an  execution  time  coefficient  of  resource  expenditure,  and  is 

a  failure  coefficient  of  resource  expenditure  for  k  =  I,  F,  and  C. 


Suppose  P^  represents  the  number  of  available  personnel,  k  =  I,  F,  or  the 
available  number  of  computer  shifts,  k  =  C.  Suppose  denotes  the  utilization 


factor  for  the  kth  resource,  [from  assumption  (i),  pj  =  1].  Then  the  effective 
available  amount  of  the  kth  resource  is  From  this  basic  formulation,  Musa52 
derives  the  following  correspondence  between  the  resources  and  the  calendar  time: 


=  MT0 


t2 

/ 

Ti 

£ 

k 


Vk 


r  ekT  +  cMki 

x  ■  =  —  dT 

L  KPK2  J 

+  Mk  '  *h)]  1 


(4.502) 


(4.503) 


where  the  index  k  can  have  the  values  C,  F,  or  I,  and  the  quantities  T^  and  T^ 


represent  the  MTBF  at  the  boundary  of  these  periods, 
values  Ti,  T2,  and  the  transition  points 


kk‘ 


V  V  ek  -  pk  p 


k  8k' 


These  boundaries  are  the 


(4.504) 


for  k,  k'  =  I,  F,  C.  The  transition  points  are  these  values  of  T  at  which  the 
derivative  of  calendar  time,  with  respect  to  execution  time  for  one  resource, 
becomes  greater  than  another.  The  resource  k  that  is  limiting  for  any  given 
MTBF,  T,  is  the  one  that  maximizes: 


V +  cpk 

pk  pkT 


(4.505) 


From  assumption  (j),  it  can  be  established  that  the  utilization  factor  for 
failure  correction  personnel  is  of  the  form: 


Pv  =  C1 


1/P 

P  *) 


F\  1/Q 


(4.506) 


where  Q  is  the  established  limitation  of  error  queue  length  (at  a  specified 
probability  P)  for  any  debugger. 


As  can  be  seen  from  the  formulation  of  the  model,  the  data  required  for 
implementation  of  the  complete  model  can  be  quite  extensive. 


Data  Requirements 
[Execution  Part] 

(a)  The  linear  execution  frequency,  f. 
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(b)  An  initial  estimate  of  the  error  exposure  ratio,  K.  (The  accuracy  of 
the  initial  estimate  can  be  low) . 

(c)  The  error  reduction  factor,  B. 

(d)  The  testing  compression  factor,  C. 

(e)  An  initial  estimate  of  the  total  number  of  errors,  N.  (The  accuracy  of 
the  estimate  can  also  be  low  since  it  is  reestimated  during  testing.) 

(f)  The  times  (measured  in  CPUs)  between  error  occurrences,  t/s. 

Execution/Calendar  Time  Part 

(g)  The  available  resources  for  both  testing  and  correction  personnel  and 
the  number  of  computer  shifts;  i.e.,  Pj,  Pj.,  and  P^. 

(h)  The  utilization  factor  for  each  of  these  resources,  i.e.,  pj  (-  1),  PF, 
and  p^. 

(i)  The  execution  time  coefficent  of  resource  expenditure  for  each  resource; 
i.e.,  0p  6F  (=0  usually)  and  0G< 

(j)  The  failure  coefficient  of  resource  expenditure  for  each  resource;  i.e., 
Pp  Hf,  and  Pg. 

(k)  The  maximum  error  queue  length,  Q,  for  a  debugger. 

(l)  The  probability,  P,  that  the  ei.ror  queue  length  is  no  larger  than  Q. 

Two  extensions  to  Musa's  Model  are  briefly  discussed  here.  The  first  appears 
in  a  paper  by  Chencweth  (Reference  57).  In  that  paper,  the  error  reduction  fac¬ 
tor,  B,  is  generalized  to  the  form  Boca^  where  Bq  is  the  initial  error  reduction 

factor  and  a  is  the  exponential  slope  of  execution  time.  Chenoweth  argues  that 
for  a  certain  class  of  software  programs,  B,  appears  to  be  exponentially  increas¬ 
ing.  The  basis  of  the  increase  is  probably  due  to  a  programmer  learning  curve 
phenomena.  The  parameter  a  can  be  estimated  from  the  relationship! 


i  n(li)  =  B0  £  ie  1 
i=l  i“l 


(4.507) 


where  n(t.)  is  the  number  of  errors  corrected  by  time  t.  for  a  specified  j  (j"l, 
i  1 

. . .  ,m  the  number  of  errors  observed)  and  Bo  is  obtained  from  a  project  of  a  sim¬ 
ilar  nature  or  using  this  relationship. 

The  second  modification  is  contained  in  a  paper  by  Musa  and  Iannino  (Refer¬ 
ence  58).  The  modification  can  actually  be  applied  to  many  of  the  previously 
considered  models,  but  it  is  illustrated  in  the  report  on  the  execution  time 
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theory  model.  The  paper  describes  a  method  of  adjusting  the  lengths  of  the  inter¬ 
vals  between  software  failures  to  compensate  for  programs  that  are  undergoing 
variations  in  length  due  to  integration  or  design  changes.  The  models  considered 
so  far  have  been  applied  to  essentially  complete  programs.  In  the  testing  pro¬ 
cess,  all  of  the  code  is  being  executed  at  one  time  or  another.  Frequently,  how¬ 
ever,  only  part  of  a  program  is  tested  and  other  parts  are  added  as  testing 
proceeds.  By  ignoring  these  variations,  estimated  MTBFs  in  the  early  stages  of  a 
project  tend  to  be  optimistic.  The  method  presented  in  this  paper  attempts  to 
account  for  the  variations  by  adjusting  the  observed  failure  intervals  to  values 
that  would  have  been  for  a  program  in  its  final  configuration  with  complete  in¬ 
spection.  The  adjusted  values  are  used  in  the  various  models  in  the  exact  manner 
as  if  they  had  been  the  actual  data.  The  reader  is  referred  to  Musa  and  lannino's 
paper  for  details. 


4.2. 12  Brooks  and  Motley’s  Models 

The  last  models  discussed  in  this  section  are  the  Binomial  Model  and  the 
Poisson  Model  formulated  by  Brooks  and  Motley  of  the  IBM  Corporation  (Reference 
59).  Their  models  try  to  account  for  the  fact  that  in  a  given  testing  period  not 
all  of  the  program  is  tested  equally,  and  in  the  development  of  a  program,  only 
some  portion  or  modules  may  be  available  for  testing.  In  addition,  in  the  cor¬ 
rection  of  discovered  errors,  additional  errors  may  be  introduced.  Each  of  the 
models  make  the  following  assumptions: 

Model  Assumptions 

(a)  The  number  of  software  errors  detected  on  each  test  occasion  is 
proportional  to  the  number  of  errors  at  risk  for  detection  which  is,  in  turn, 
proportional  to  the  remaining  number  of  errors. 

(b)  The  proportionality  factor  or  probability  (denoted  as  q  for  the 
binomial  model,  and  4>  for  the  Poisson)  of  detecting  any  error  during  a  speci¬ 
fied  unit  interval  of  testing  is  constant  over  all  occasions  and  independent 
of  error  detections. 

(c)  The  errors  reintroduced  in  the  correction  process  are  proportional 
to  the  number  of  errors  detected. 

For  their  formulation,  Brooks  and  Motley,  develop  the  models  both  for  a  module 
application,  in  case  only  module  testing  is  done,  and  for  the  entire  program  sys¬ 
tem  testing. 


4.2.12.1  Binomial  Model  (module).  Suppose  a  module,  the  jth,  from  the  pro¬ 
gram  is  given  for  testing  for  the  first  time.  Then  the  expected  number  of  error 
occurrences  in  that  module  in  the  first  unit  interval  of  the  test  occasion  is: 

n^  =  WjNq.  (4.508) 
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This  is  obtained  from  assumptions  (a)  and  (b), 


where  w.  is  the  weight  assigned 
3 


to  module  j ,  N  is  the  total  number  of  errors  in  the  system  at  the  beginning  of 
testing,  and  q  is  the  error  detection  probability  given  in  assumption  (b) .  Brooks 
and  Motley  define  a  test  occasion  as: 


"an  event  of  error  data  collection;  each  occasion  should  have 
a  time  interval  associated  with  it;  otherwise,  the  implica¬ 
tion  to  the  model  is  that  all  test  occasions  are  of  equal 
length  of  time....  One  additional  important  assumption  made 
here  is  that  one  occasion  be  comparable  to  every  other  occa¬ 
sion  in  terms  of  the  time  spent  (testing  effort)  in  detecting 
errors." 


The  weight  factor  can  be  taken  as  the  ratio  of  the  size  of  the  module  (as  measured 
by  number  of  lines  of  source  code  or  object  program  size)  to  the  total  program 
size. 


For  the  second  unit  interval  of  testing  on  the  jth  module,  the  expected  num¬ 
ber  of  errors  to  be  detected  is: 


[WjN  -  WjNqlq  =  [WjN(l  -  q)]q.  (4.509) 

There  were  w.Nq  expected  errors  in  the  first  unit  interval  of  time  leaving  w.N  - 

J 

w.Nq  errors  subject  to  detection  in  the  second.  Thus  the  expected  number  of 

J 

errors  in  the  second  time  interval  is:  [number  of  errors  subject  to  detection] 

*q  =  [w.N(l  -  q)]q.  In  general,  for  the  ith  unit  interval  of  time  in  the  first 
J 

testing  session,  the  expected  number  of  errors  is: 

w.N(l  -  q)i-1q.  (4.510) 

J 

The  total  number  of  errors  expected  for  the  entire  first  testing  occasion  is  then: 


w.Nq 
J  H 


U 


(4.511) 

(4.512) 

(4.513) 


where  is  the  number  of  unit  test  intervals  making  up  the  first  test  occasion, 

or  the  total  test  effort  expended  on  module  j  during  the  first  test  occasion,  and 
K-  . 

„  -  fi  »  ri  _  i 
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When  module  j  is  tested  for  the  second  test  occasion,  the  number  of  errors  at  risk 
in  module  j  is: 

WjN  -  n^  +  rn^  (4.514) 

where  n^  is  the  number  of  errors  detected  in  the  first  testing  period  and  rn^  is 

the  number  introduced  into  the  program  as  a  result  of  correcting  those  n^  errors, 

(assumption  (c)).  The  total  expected  number  of  errors  in  the  second  test  period 
can  be  shown  as  was  done  for  the  first  to  be; 


n2j  =  -  “ij^j 

where  a  =  1  -  r,  (the  probability  of  correcting  code  without  introducing  new 
errors)  and 


q2j  =  (1  -  (1  -  q)  2j )  . 


(4.515) 

In  general  for  the  ith  testing  period,  the  expected  number  of  errors  detected  is: 

n.  .  =  (w.N  -  ON.  ,  .)q..  (4.516) 

ij  J  i-l. J  HiJ 

(4.517) 


where 


=  N. .q. . 


i-1 


N.  .  .  =  Y'  n  . 

1-1  m=i 


(4.518) 


(the  total  number  of  errors  found  up  to  the  ith  testing  period), 


K. . 


-  (i  -  (i  -  .)  1J)  , 

(4.519) 

and 

N.  .  =  (w.N  -  <»N.  ,  .)  , 
iJ  J  i-l.  J 

number  of  errors  remaining  in  the  jth  module. 

(4.520) 

the 

One  notices  that  IL  j ,  the 

amount  of  testing  effort  expended  in  the  ith  period,  can  be  different  from  one 
testing  period  to  the  next.  The  only  restriction  is  that  the  probability,  q,  of 
detection  for  any  error  is  the  same  from  period  to  period.  This  means  that  the 
times  can  vary  for  each  testing  period,  but  the  testing  approach  should  be  the 
same. 


Brooks  and  Motley  establish  the  MLEs  of  the  three  unknowns,  (N,  q,  a)  of 
their  model  as  the  solutions  of  the  following  equations. 
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Estimates  -  Maximum  Likelihood 


K  J  r 

=  E  £  k; 

i=l  j=l L 


"id  ■  “id 


8&nL  _  n  _ 

aq  " 


K  j  r  n..K.. 

E  E,  / - 3  3  kT." 

1-1  l3( 


+  w.K. .£n(l  -  q) 
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(4.521) 


(4.522) 
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M 

where  the  likelihood  function 


.£n(l  -  q)  (4.523) 


K  J  / N..\  n..  N..  -  n.. 

=  n  n  (  XJ  J  q,.  3  (i  -  q)  J  J  > 

i=i  d=i  \»(H/  1J 


(4.524) 


K  =  the  number  of  test  occasions,  J  =  number  of  modules  in  the  system, 
and  n.  .  is  the  actual  number  of  errors  observed  on  the  ith  testing  occasion  of  the 

jth  module.  These  equations  do  not  have  a  solution  if  -  n^  becomes  negative. 

It  could  happen  that  the  effective  number  of  errors  at  risk,  Nij ,  becomes  smaller 

than  the  actual  number  of  errors  observed  for  the  jth  module  on  the  ith  testing 
occasion.  In  that  situation,  it  is  recommended  that  the  system  model  be  applied. 

4.2.12.2  Poisson  (module).  As  in  the  Binomial  Model,  suppose  =  (wjNj  “ 

aN.  *  .)  is  the  effective  number  of  errors  at  risk  in  module  j  at  the  beginning 
J  )  J  ^ 

of  the  ith  testing  period.  Using  assumption  (b),  the  expected  error  detection 
rate  for  the  first  unit  interval  of  length  t  is  N.  .<{>.  Thus  the  expected  number 
of  errors  that  are  detected  is  the  error  detection  rate,  N^,  times  the  length  of 
the  testing  interval;  t,  i.e., 

“id  =  V-  (4-525) 

At  the  end  of  the  first  unit  time  interval  testing  period  of  module  j,  the  number 
of  errors  remaining  is: 

N.  .  -  N.  .(|>t  =  number  of  errors  in  the  module  at  the  beginning  of  the 
XJ  first  period  minus  the  number  of  errors  detected  during 


the  first  testing  period. 


(4.526) 
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The  error  detection  rate  for  the  second  testing  period  is: 

[N.j  -  N.^tU  =  {N.jO  -  <>t)]4»  (4.527) 

so  the  expected  number  of  errors  is  therefore: 

[N±J(1  -  <j>t)]$t 

for  the  second  unit  interval  of  testing.  In  general,  for  the  t..  unit  interval  of 
testing  of  the  jth  module  in  the  ith  testing  period:  1‘* 

(a)  The  number  of  errors  remaining  at  the  beginning  of  the  interval 

t.  .-1 


=  N.  .(1  -  40 

ij  T 


ij 


(b)  The  expected  error  detection  rate  is: 

t.  .-1 

=  N..(l  -  40  1J  4 

and  thus, 

(c)  The  expected  number  of  errors  detected  is: 


(4.528) 


(4.529) 


N^d  -*40 


t.  .-1 

<J>t  . 


(4.530) 


The  length  of  the  unit  testing  interval,  t,  is  then  normalized  to  1  (example:  1 
day,  1  week,  etc).  Thus  the  total  number  of  expected  errors  for  the  jth  module 
on  the  ith  testing  occasion  is: 


t.  . 

ij  M 

T 

n.  .  = 

^  2=1 


N.^d  -  4) 


2-1 


where 


=  N.  .4.  • 
ij  ij 


4Aj  =  1  -  (1  -  4) 


t. . 

ij 


The  likelihood  function  is  then: 

n. .  -  N. .4. . 

K  J  (N..4-.)1Je  ^ ^ 

L  =  n  n 

i=i  j=i 


n.  . ! 

ij 


The  MLEs  are  then  obtained  as  the  solution  to  the  equations: 


fiiii 

(4.531)  ' 

* 

**  * 

tt*  ~-2rj 

(4.532) 

V*V 

-  k  •, 

(4.533) 
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Estimates  -  Maximum  Likelihood 


9£nL 

9N 


K  J  Tn..  “I 

=  0  =  £  £  wj  iH-  “  ^ij  I 
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=  o=  £  £  t  (l  -  4.)  1J  - 1J  t  -  -  N 

i=l  j=l  J  L1  '  JJ 


(4.535) 


(4.536) 
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(4.537) 


4.2.12.3  Binomial  (system).  For  this  model,  the  overall  program  is  con¬ 
sidered  as  a  whole.  Keeping  the  same  notation  as  was  used  in  the  Binomial  (mod¬ 
ule)  paragraph,  if  JL  is  the  index  set  of  those  modules  tested  on  occasion  i , 

then  the  total  number  of  errors  remaining  in  the  program  and  subject  to  detection 


si  =  £  • 


(4.538) 


Since  the  system  is  being  considered  as  a  whole,  the  test  effort  involved  for  the 
system  can  be  considered  as  a  whole  rather  than  on  a  modular  basis.  The 


q£j  =  [1  -  (1  -  q)  1J] 
of  the  modular  section  is  then  replaced  by: 


(4.539) 


q.  =  [1  -  (1  -  q)  *] 


(4.540) 


where  is  the  system  test  effort  (e.g.,  computer  CPU  time)  expended  ou  the  ith 

test  occasion.  Combining  the  previous  information,  the  total  number  of  expected 
errors  in  the  system  for  the  ith  testing  occasion  is: 


n.  =  N.q,  . 
l  in 

The  likelihood  equation  for  the  system  is  therefore: 


K  /N.  \  n.  N.-n. 

L  =  n  1  q*  a  -  q.)1  1 

1=1  w 


(4.541) 


(4.542) 


i  *f  *  •r\  ''X  # . »  * . 

k.V.V.A-V.V 


where 


j 

n.  =  53  n.  .  =  total  number  of  errors  found  in  the  program 
1  j=l  on  the  ith  testing  occasion. 


(4.543) 


The  MLEs  are  then  obtained  as  the  solutions  to  the  following  set  of  equations: 
Estimates  -  Maximum  Likelihood 


8£nL 

3N 


=  0  = 


K.£n(l  -  q)\  £  w. 
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(4.545) 
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+  K.An(l  -  q) 


(4.546) 


Again  if  IT  -  n.^,  the  difference  between  the  expected  number  of  remaining  errors 

and  the  actual  number  of  errors  found  on  the  ith  testing  occasion  is  negative,  no 
solution  exists  to  these  equations. 

4.2.12.4  Poisson  (system).  Using  the  expression  for  the  number  of  errors 
at  risk  at  the  start  of  the  ith  testing  occasion  for  the  system,  i.e., 


S.  =  L  (w.S-ON.  .  .)  , 

1  jeJ.  J  i-'-J 

the  to^al  expected  number  of  errors  is: 


where 


n.  =  N.<j>. 

X  11 


t. 

4>.  =  l  -  (l  -  <|»)  1 


(4.547) 


(4.548) 


(4.549) 


and  t^  is  the  total  time  spent  for  the  ith  testing  occasion. 
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The  likelihood  function  is  therefore: 


L 


K 

=  n 

i=l 


ni  'Ni^i 
(N .  0 . )  e 


iTi* 


n. ! 
1 


(4.550) 


where  n.  is  the  total  number  of  errors  detected  during  the  ith  testing  occasion. 
i 

The  MLEs  are  again  obtained  as  the  solutions  to  the  following  system  of  equations: 


Estimates  -  Maximum  Likelihood 


and 


3£nL 

3N 


V1 


(1  -  4)fci 


(4.551) 

(4.553) 


(4.554) 


The  various  sets  of  equations  given  in  the  previous  paragraphs  can  all  be 
solved  using  the  Newton-Raphson  method  with  the  warning  ab_  t  the  (JiL  -  n^'s  be¬ 
coming  negative  applying  as  in  the  binominal  formulations. 


All  of  the  models  by  Brooks  and  Motley  were  applie  ■  to  real  life  and  simu¬ 
lated  data.  One  criticism  that  might  be  made  against  their  models  is  the  assump¬ 
tion  of  a  constant  detection  probability,  q.  In  a  testing  environment,  the  usual 
situation  is  that  q  varies  over  time.  This  is  due  to  limiting  resources,  the 
easier  errors  are  found  at  the  beginning,  while  the  hidden  errors  are  discovered 
much  later  and  at  greater  effort,  and  there  is  a  learning  curve  effect  on  the 
testers.  Brooks  and  Motley  do  however  try  to  account  for  this  by  considering  an 
extension  to  their  basic  models.  They  allow  for  the  probability  of  detection  to 
increase  at  a  constant  amount  until  it  reaches  a  point  where  it  levels  off .  The 
resulting  equations  for  this  extension  are  very  complex  and  would  be  difficult 
to  implement  on  a  computer.  The  reader  is  referred  to  their  paper  for  additional 
details . 


4.3  BAYESIAN  MODELS 

The  class  of  models  considered  in  this  paragraph  formulated  software 
reliability  modeling  in  a  Bayesian  framework.  The  models  employ  a  "subjective 
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approach  to  the  meaning  of  software  reliability  in  constrast  to  the  traditional 
"frequentist"  approach.  Previous  models  only  allowed  for  change  in  the  reli¬ 
ability  of  a  program  whenever  an  error  was  discovered  and  subequently  corrected. 
Bayesian  models  take  the  subjective  viewpoint  that  as  the  software  is  tested,  if 
no  errors  are  discovered,  there  is  more  confidence  in  the  program  and  this  is 
reflected  in  increasing  reliability.  The  reliability  of  a  program  should  be  a 
reflection  of  the  number  of  errors  discovered  and  the  length  of  error-free  testing 
time  periods. 

Another  important  argument  given  in  support  of  i  Bayesian  approach  deals  with 
counting  errors.  All  of  the  models  considered  so  tai.  assume  that  the  hazard  rate 
function  is  directly  proportional  to  the  number  of  errors  in  the  program  at  the 
time.  From  this  assumption,  it  is  directly  determined  that  the  reliability  is  a 
function  of  this  count.  This  is  the  reason  the  models  considered  in  the  previous 
paragraph  are  concerned  with  estimating  this  total.  The  Bayesian  approach  argues 
that  a  program  with  two  or  more  errors  in  little  exercised  portions  of  code  is 
considered  more  reliable  than  one  with  only  one  error  in  a  frequently  executed 
section  of  code.  The  estimation  of  the  total  number  of  errors  present  can  be  of 
use  to  the  software  manager  in  making  determinations  of  resource  allocation,  but 
it  should  not  be  the  driving  factor  in  reliability  considerations.  One  should  be 
concerned  with  measuring  operational  reliability. 

A  number  of  models  which  attempt  to  do  this  are  now  considered. 


4.3.1  Littlewood's  Bayesian  Debugging  Model 

The  first  model  considered  within  this  class  was  proposed  by  Bev  Littlewood 
of  the  City  University  of  London.60’61’62’63’64  The  model  reformulates  the 
Jelinski-Moranda  Model  (Paragraph  4.2.3)  into  a  Bayesian  framework.  The  Jelinski- 
Moranda  Model  postulates  that,  at  any  point  in  time,  the  error  rate  is  propor¬ 
tional  to  the  number  of  errors  remaining  in  the  program.  This  is  expressed  as  for 
any  time  t,  for  t^_j  <  t  <  t^, 

Z(t)  =  <|>(N  -  i  +  1)  (4.555) 

where  the  t^'s  are  the  times  of  error  occurrences.  By  making  the  assumption  that 
the  times  between  error  occurrences,  i.e.,  X^  =  t^  -  t^_j,  follow  an  exponential 
distribution,  the  probability  density  function  for  is  seen  to  be: 

f(X.)  =  <|>(N  -  i  +  l)expl-<j)(N  -  i  +  1)X±] .  (4.556) 

The  model  inherently  makes  the  assumption  that  all  errors  contribute  equally; 
namely  <|>,  to  the  overall  error  rate.  The  Bayesian  viewpoint  objects  to  this 
assumption.65  Each  error  does  not  contribute  equally  since  the  correction  of 
errors  in  the  beginning  of  the  testing  phase,  has  more  of  an  effect  on  the  program 
than  ones  corrected  later.  Again  the  argument  that  a  program  with  two  errors  in 
rarely  exercised  code  is  more  reliable  than  a  program  with  only  one  error  in’  a 
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frequently  exercised  section  surfaces.  All  errors,  therefore,  do  not  contribute 
equally.  Littlewood  postulates  that  the  error  rate 

=  Z(t)  =  4>(N  -  i  +  1),  t._1  <  t  <  (4.557) 

should  be  treated  as  a  random  variable,  not  as  a  constant.  By  assuming  that  the 
remaining  errors  have  different  occurence  rates  0i,  4*2 »  the  overall  fail- 

ure  rate  is  then: 

\  =  <t»l  +  <h  +  •••  +  ^N-i+l  *  (4.558) 

By  treating  the  <j)^'s  as  random  variables  (since  it  is  not  known  what  they  are), 
the  overall  rate  as  a  random  variable  is  obtained.  (Notice  that  if  all  of  the 
<|k's  are  assumed  to  have  a  degenerate  distribution  at  the  point  <{>,  i.e.,  <f>i=<)>2 
=  ...  =  =  $  with  probability  1,  then  =  <|>(N  -  i  +  1)).  The  specific  as¬ 

sumptions  for  this  model  are: 


Model  Assumptions 


(a)  The  individual  failure  rates  of  the  errors  in  the  program  are  assumed  to 
be  independent  random  variables  each  with  a  prior  distribution  that  is  assumed 
gamma  with  parameters  a  and  p,  i.e., 


gC^)  =  g(<h)  = 


eV1  e~f 


,  4>  >  o 


(4.559) 


for  all  i  and  j . 

(b)  For  a  given  error  rate,  the  time  between  error  occurrence  X^  =  t^ 


ti_^  is  assumed  to  be  exponential  with  mean  1/X^ ;  i.e., 


f(X.  \.)  =  X.e 

li  l 


-X.X. 

l  i. 


X.  >  0  . 

,1 


(c)  X  =  <h  +  <{>2  +  . . .  +  <t>N_i+1 


(4.560) 

(4.561) 


after  i  -  1  errors  have  been  detected  and  corrected. 


(d)  When  a  software  error  is  detected,  it  is  immediately  corrected  without 
the  introduction  of  additional  errors. 

(e)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

The  model  is  developed  as  follows.  Suppose  Xi,...,  X  are  the  times  between 
errors  occurrences,  i.e.,  ‘ 


X.  =  t.  -  t.  . 
i  i  l-l 
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(Preferably  is  measured  in  CPU  time  rather  than  wall  clock  time.)  At  the  time 
the  ith  error  is  discovered  and  corrected  from  assumption  (c) , 

\+1  =  4>i  +  4*2  +  •••  +  <1^  •  (4.563) 

Now  suppose  the  occurrence  rate  ^  for  any  one  of  the  remaining  N  -  i  errors  is 
considered.  The  density  function  for  <j>k  is  pdf (^^  |  given  that  error  was  not  found 
in  (0,t)),  where  t  is  the  current  testing  time 

P{no  failure  by  that  error  in  (0 , t) |  $k  =  4>k}  pdf((|)^) 


=  n  _  i  *>  ^ _  tv 

“  J  P{no  failure  by  that  error  in(0,t)  <J>k  =  (f)^}  pdf(<j»^)d^ 


(4.564) 
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(4.565) 


(P+t)0  <J.kC(“1e“(P+ti)<t,k 

TO) 


(4.566) 


Thus  <]>k  has  a  density  function  that  is  also  gamma  with  parameters  a  and  p+t. 
Since  is  a  sum  of  independent,  identically  distributed  random  variables,  \^+1 
is  also  gamma  with  parameters  (N  -  i)a  and  P+t.  Thus,  the  unconditional  distri¬ 


bution  of  the  time  to  the  next  failure  X.,,  is: 

l+l 

f(xi+i)  »  /*  £«1+i|  W 


(4.567) 
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_  (B  +  t)(N“i)a  F  f (N  -  i)ct  +  1] 

'  HTN  -  i)af  *  ^  +  ;  +  Xi+i)(H-i)'«+l  ‘ 


(4.570) 
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(4.571) 


This  is  a  Pareto  distribution.  From  this  basic  result,  Littlewood  derives  a  num¬ 
ber  of  quantities.  The  reliability  function  after  i  errors  are  discovered  is 
found  as: 


R(x)  =  jXi+1  >  x) 

*  1  '  PtXi+1  <  x) 

X 

■ 1 '  /  £(WdXiti 

_[  (P  +  t)  1  (N  -  i)a 

[(P  +  t  +  x)J 


The  failure  rate  function  is  then  obtained  as: 


Z(x)  = 


R'(x)  _  (N  -  i)q 
TtxT  ~  (P  +  t  +  x)  * 


(4.572) 

(4.573) 

(4.574) 

(4.575) 

(4.576) 


Thus,  the  failure  rate,  immediately  after  i  errors  have  been  discovered  and  t 
amount  of  testing  has  been  employed,  is: 


Z(0)  = 

P  +  t 


(4.577) 
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Notice  how  this  unconditional  failure  rate  changes  after  testing  progresses.  As 
t  gets  larger,  the  hazard  rate  decreases  reflecting  the  increased  confidence  in 
the  program.  The  hazard  rate  also  decreases  whenever  an  error  is  discovered  and 
corrected. 

The  MTBF  is  found  from  the  Pareto  distribution  as: 

MTBF  =  E(X1+1)  =  /  X1+1f(X1+1)dXi+1  (4.578) 
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8  +  t 

(N  -  i)a  -  1 


(4.579) 


which  exists  as  long  as  (N  -  i)a  >  1. 

From  this  basic  model  formulation,  Littlewood  approaches  the  problem  of  pre¬ 
diction  of  future  reliability  in  one  of  two  ways.  The  reliability  of  the  program 
can  be  estimated  after  some  specified  execution  time  has  elapsed  or  after  some 
specified  number  of  errors  has  been  removed.  The  development  of  these  two  ap¬ 
proaches  is  provided  in  Reference  61  and  is  not  repeated  here;  however,  the  uses 
of  those  results  are. 

For  the  first  approach,  suppose  t^  amount  of  testing  is  performed  and  i 

errors  are  discovered  and  corrected.  Now  suppose  At  additional  amount  of  testing 
is  done.  Then  the  reliability  of  the  program  at  the  time  t^  +  At  is  shown  to  be: 


R<*> "  [l  ) 


P  +  t. 

P  +  +  At  +  x 


(4.580) 


From  this  relationship,  the  amount  of  additional  testing  needed  in  order  to 
achieve  a  target  reliability  can  be  determined.  If  the  desired  reliability  is  r 
for  a  specified  error-free  run  time  of  xo,  then  the  additional  testing  time  re¬ 
quired  is  the  value  At  that  solves  the  equation: 


=[  1  -(e  +  A t  j 


P  +  t. 

p  +  t^  At  +  xq 


otl  N-i 


(4.581) 


Littlewood  also  shows  for  this  approach  that  the  required  additional  testing  time 
At  to  achieve  a  specified  target  failure  rate,  \q>  is: 


for  -  iJaCMtJ0 


‘  (P  +  t.) 


(4.582) 


For  the  second  approach,  suppose  i  errors  are  observed  and  corrected.  In¬ 
terest  lies  in  the  times  between  error  occurrence  of  the  next  k  errors,  i.e., 
xi+j »  j*l,...,k.  Littlewood  first  derives  the  distribution  of: 


" A(P  +  t.) 

=  (N  -  i  -  k)« 


(4.583) 


where  A  is  the  failure  rate  at  the  occurrence  of  the  i  +  k  error.  From  the  pre¬ 
vious  results,  the  following  is  obtained: 


(N-i  -  k)q 
i+k 

p  +  t  +  Ex. 

j=i+l 


(4.584) 
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The  distribution  of  Z  is  beta  with  parameters  N-i-k+1  and  k,  so  that  the  expected 
value  of  A  can  then  be  obtained  as: 


ElA>  ’  (”  'e  i  t. 


k)a  (N  -  i  -  k  +  l)a 

(N  -  i  -  k  +  l)o  +  1 


(N  -  i) 


(N  -  i)d  +  1 


(4.585) 


From  this  result,  the  number  of  additional  error  corrections,  ko,  that  are  re¬ 
quired  to  make  the  E{A}  -less  than  a  desired  level  Ao  can  be  established.  This 
is  the  smallest  integer  ko  satisfying: 

(N-  ~ .  \  Z  .  (N_-_i  -  k  +  l)q _  .  .  .  (N  -  i)__  <  x 

p  +  tt  (N  -  i  -  k0  +  l)a  +  1  (N  -  i)o  +  1  A°  * 

(4.586) 

It  might  also  be  asked  how  many  additional  error  corrections  are  necessary  in 
order  to  be  at  least  Y  percent  certain  that  A  <  Ao*  This  is  the  smallest  integer 
ko  such  that: 


Z  < 


r  o  + 1. 

[W-'  i  - 


) 

k0)a 


-\ot 


>  Y 


(4.587) 


where  Z  is  from  a  beta  distribution  with  parameters  N-i+ko+1  and  ko 


For  this  model,  there  are  three  unknowns:  N,  a,  and  p.  They  can  be  estimat¬ 
ed  using  the  maximum  likelihood  procedure  or  least  squares.  If  X^,  i  =  l,...,n  is 

the  time  between  error  occurrences,  then  from  the  assumptions,  the  likelihood 
function  is: 


n 


k(N,oi,p)  -  n  f(X.  x.  . x.) 

i=l  1  1  1 


n 


iUj  (N  -  i  +  l)a(p  +  t^j)^  "  1  +  1)0 
(P  +  Vl  +  X.)^“  1  + 


(4.588) 


(4.589) 


where 


‘‘-1  =  &  "i 

is  the  time  of  occurrence  of  the  (i  -  l)st  error. 

A  A  A 

The  MLEs  N7  ,  o\.  ,  and  p^  are  the  estimates  which  maximize: 


(4.590) 


jg 

Mi 

ty, 

fcv 
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Estimates  -  Maximum  Likelihood 


L(NL,0L,PL)  =  max  L(N,a,p)  . 
N,a,p 


(A. 591) 


Littlewood  points  out  that  this  maximization  search  can  be  restricted  to  the  two- 
dimensional  space  of  N  and  p  as: 


aL  = 


'  A 

n-i  rsT  +  t.  "1  -  r pT  +  t  ' 

£  jd  3 - I  L  k  £n  X— S 

i=l  Lftr+tJ  L  L  Pr  J 


*  *  *-T.  nJ  *-  T,  J 

The  least  squares  estimates  are  those  N,  a,  and  p  which  minimize 

S(N,  L  ‘ 

L  A  A  A 

using  equation  (4.579).  The  least  squares  estimates  NT  TO>  pT  TC,  and  aT 
chosen  so  that:  L,Lw  L,LS  L’ 

S(NL,LS*  PL,LS»  “l.LS5  =  (^p)’  **  P)  ( 

and  are  found  as  the  solution  to  the  following  system  of  equations: 

9N  "  °»  8a  "  °»  and  9p  "  °*  ( 

The  least  squares  estimates  are  then  the  solutions  to  the  equations: 
Estimates  -Least  Squares 


t  xi  =  f  *  Vi 

1=1  (ni,ls  •  1  +  "Xis  -  1  ’  1=1  1(nl,ls  -  1  +  »%is  -  112 


Xi«L,LS  +  ‘l-l> 


n 

-  z 


^L.LS  +  tl- 1  >* 


&  l»L.LS  -  1  +  1)aI,LS  ’  *!2  1=1  ICl.LS  -  1  +  1)0l,LS  ‘  11 


(4.596) 


(4.597) 


(4.592) 

■VO, 

Pk  *.*•- 

(4.593) 

[vVj 

,LS  are 

$ 

m 

(4.594) 

-  * 
,v 

v*| 

K*-  *7) 

(4.595) 
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and 

A  VVlS  -  *  +  15  ‘VlS  +  Vl>  £  (iiI,LS  -  1  +  1)(Vl,S  *  VP’ 

£  ioC  - 1  + 1)  C  -  u‘  ■  =  ^ 


i=1  ^(NL,LS  “  1  +  1)aL,LS  "  1^3 


(4.598) 


To  implement  this  model,  the  data  required  are: 


Data  Requirement 


The  times  between  error  occurrences,  i.e.,  the  XJs,  or  the  times  of  error 
occurrences,  i.e.,  the  t/s  where 


X. 


l 


t. 

l 


Vr 


(4.599) 


Once  the  parameters  N,  a,  and  p  are  estimated,  all  of  the  previous  quantities 
developed  in  this  paragraph  can  be  estimated  by  replacing  the  parameters  with 
their  corresponding  estimates. 


An  alternate  Bayesian  modification  of  the  Jelinski-Moranda  Model  is  given  in 
a  paper  .by  Littlewood  and  Sofer.66  In  that  paper,  the  times  between  error  occur¬ 
rences,  i.e.,  X^’s,  are  assumed  exponential,  but  with  parameters 

A.  s  X  -  (i  -  1)4>  i  =  l,.v.  ,n.  (4.600) 

Constrast  this  with  the  model  considered  in  this  paragraph  of 


Ai  *  4>l  +  $2  +  *  •  +  ^-i+1  (4.601) 

For  both  formulations,  the  X.'s  are  taken  as  random  variables.  For  the  alternate 
model,  the  X  and  are  taken  as  independent  random  variables  with  prior  distri¬ 
butions  Gamma(b,c)  and  Gamma(f,g),  respectively.  All  of  the  quantities  devel¬ 
oped  in  this  section  for  the  first  Bayesian  Model  are  developed  for  the  analogous 
ones  in  Littlewood  and  Sofer’ s  report.  They  are  not  repeated  here. 


4.3.2  Littlewood  and  Verrall's  Bayesian  Reliability  Growth  Model 

The  next  model  considered  is  the  Bayesian  Reliability  Growth  Model  proposed 
by  Littlewood  and  Verrall. 67,68,69  The  model  tries  to  account  for  error  genera¬ 
tion  in  the  corrective  process  by  allowing  for  the  probability  that  the  program 
could  be  worsened  by  correcting  the  error.  The  intention  is  to  make  a  program 
more  "reliable"  when  an  error  is  discovered  and  corrected,  but  there  is  no 
assurance  that  this  goal  is  achieved.  With  each  error  correction,  a  sequence  of 
programs  is  actually  generated.  Each  is  obtained  from  its  predecessor  by 
attempting  to  correct  an  error.  Because  of  the  uncertainty  involved  in  this 
correction  process,  the  relationship  that  one  program  has  with  its  predecessor 
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cannot  be  determined  with  certainty.  This  is  a  second  source  of  uncertainty  in 
the  modeling  of  software  reliability  (the  first  dealing  with  the  variation  of  the 
input  to  the  program).  The  specific  assumptions  for  the  model  are: 


Model  Assumptions 

(a)  Successive  execution  times  between  failures,  i.e.,  X.,  i=l,...,n,  are 
independent  random  variables  with  probability  density  functions  1 

-X.X. 

fCXilV  =  A^e  1  1  >  D  .  (4.602) 

That  is  X^  is  assumed  exponential  with  parameter  X^. 

(b)  The  X^'s  form  a  sequence  of  independent  random  variables  each  with  a 
gamma  distribution  of  parameters  a  and  t|>(i),  i.e., 

a  a-1  -i{f(i)X. 

[4i(i)]X.  e 

g(\.)  =  - - -  \.  >  0  .  (4.603) 

1  r(a)  1 

The  function  tj<(i)  is  taken  to  be  an  increasing  function  of  i  that  describes  the 
"quality"  of  the  programmer  and  the  "difficulty"  of  the  programming  task.  A  good 
programmer  should  have  a  more  rapidly  increasing  function  4*  than  a  poorer  pro¬ 
grammer.  The  .function  reflects  past  and  future  changes  in  reliability  as  a 
growth  process. 

(c)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

By  requiring  the  function  to  be  increasing,  the  condition 

P{X(j)  <  A}  >  P{A(j  -  1)  <  £}  (4.604) 

for  all  j  is  satisifed.  This  reflects  that  it  is  the  intention  to  make  the  pro¬ 
gram  better  after  an  error  occurs  and  is  corrected,  but  it  cannot  be  assured  that 
our  goal  is  achieved. 

When  the  two  sources  of  randomness  are  put  together,  then 


f lxii a ,  4>(i)l  = 


/  f(xi|V  8(Xi)dXi 
/  - 


-Vi 


a  ,a-l 


e  xdX. 

i 


(4.605) 


(4.606) 


.'■j'S 

£ i 

An 

£i 
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of  [V(i)]q 

tx±  +  *jj(i))a  +  1 


x.  >  0. 
x 


(4.607) 


Notice  that  the  x.'s  are  no  longer  exponential.  They  now  have  the  Pareto  distri¬ 
bution.  The  joint  density  for  the  x/s  is  then 

a11  n  [«|i(i)]a 

£[xi . xnl°-  *(t)]  =  - - 1:1 -  *!  >  0.  i  =  1 . n  . 


n  [x  +  Ui))a+1 

i=l  1 


Littlewood  and  Verrall  suggest  the  following  forms  for  the  t|f  function: 

4i(i)  =  po  +  Pli  (linear) 

and 

i|i(i)  =  Po  +  Pli2  •  (quadratic) 


(4.608) 

(4.609) 

(4.610) 


(Littlewood  finds,  for  one  set  of  data  on  which  the  model  was  applied,  that  the 
linear  function  is  superior  to  the  quadratic  function.68)  In  either  case,  the 
likelihood  function  is  now  a  function  of  three  unknowns  (a,  po,  and  pi).  MLEs 

A  A 

could  be  found  by  finding  the  a^y,  Pq  ^y,  and  pj  ^y,  which: 

A  A  A 

L  P0  lv*  LV^  =  max  L^a’  P°»  M  (4.611) 

*  *  («»  Po.  Pi) 


where 


L(a,  Po,  Pi)  =  f(xi,...,xnj  a,  Po,  Pi)  . 

These  MLEs  are  the  solutions  to  the  following  system  of  equations: 


(4.612) 


v.V.-yv--./-. 

w  w 

KSH  t  ' 

«  £3 


kvc 

IVW! 
»v  Wi 

W.vj 

V-:v 


P  fa 

D'S*; 

K\'>! 

i*.  J  ■**  A  *»  , 

L«l  \/ 

(Vt 


i.i 

,V  r 

a 

r.;i  L'','/ 

k,  •!  .  7 

y.j  kV... 

1  I.V.* 

LV;*: 


H  E 

B  m 


:  *  vf 

t.'1,  h"’*V 

ft  fly.-' 


Estimates  -  Maximum  Likelihood 

^  ^  n 

oL  _  n  .^p  \  V' 


g-J  +  E  to  *(i)  -  L  An  (x.  +  iKi)]  =  0 


ai  _ 


-a  E 


i=l  *(i) 


(a  +  1)  E 


i=l  xt  +  ijf(i) 


(4.613) 

o  ■■  r»  ; 

L'*" 

f'*  I-1.',*.’ 

L-, 

k‘.  **, 

•  »  •  •  % 

(4.614) 

\s  \V  s 

i  Li 


1’-". 

EL. .... 
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ar  n  n  •  ' 

1  =  8  E  -  (S  +  1)  £  - — 

P*  i=l  t|j(i)  i=l  +  4*(i) 


4»(i)  =  Po  +  Pii  or  p0  +  Pli2 


(4.615) 


(4.616) 


o  y, 

jj  - 


i'  =  i  or  i2  . 


(4.617) 


Littlewood  and  Verrall  eliminate  the  parameter  a  through  a  Bayesian  analysis.  By 
assuming  a  uniform  prior  for  a,  it  can  be  shown  (Reference  67)  that  the  distribu¬ 
tion  of  x.  is: 
i 


The  MLEs  for  Po  and  Pi  are  those  parameters  which  then: 


(4.618) 


i  k 

(:  C} 

i  d 


5: 

H  K> 

\\ 

V-# 


ty 


UPo.Pi)  =  max  L(Po>PO  s  max  n  f(x,  p0,Pi)  • 
(Po.Pi)  (Po.Pi)  i*l 


(4.619) 


Littlewood  and  Verrall  present  an  alternative  way  of  estimating  Po  and  Pi  based 
upon  goodness-of-fit .  The  reader  is  referred  to  their  paper67  for  details. 


Another  procedure  for  estimation  is  based  upon  least  squares.  Since 


fixja,  *(i) ]  =  — gutfiyi— 

11  [X,  +  4(i)l  1 


the  MTBF  is 


JP  a  x.  [iKi)l 

E{X  }  =  /  - i - 5CT  dx  , 

1  *  lx.  +  4(01 

=  i mi 

a  -  1 


(4.620) 


(4.621) 


(4.622) 


i«\  4*. 

!.\  v' 


.  i'\ 

i  e 


provided  a  >  1. 

The  least  squares  estimates  are  those  parameters  which  minimize: 


s<«,  Po.  Pi)  =  ±L  -  I*£i21  V 

i=l  \  1  a  -  1  / 
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,„«•  vW  OUSL*.  J0*  -"Ui  iMUliii  «JUk  ,.i 

$ 

3 

iS 


In  particular,  if 


<J>(i)  =  Po  +  Pi* 


(4.624) 


the  least  squares  estimates  aLg,  pQ  Lg,  and  ^  Lg  satisfy  the  following  system  of 
equations :  * 


(4.625) 


(4.626) 


(4.627) 


Estimates  -  Least  Squares 

n  a 

I>2  Ci) 

H  =  Sis  '  1  '  - -  ' 

E*i  ici) 

i=l 

A  A 

as  A  x  11  Pq.ls  _  pi.Lsn(n  *  ])  _  0 
1-1  Xl  '  «ls  -  l'  2<aIg  -  1) 


n 

3S_  _  r  (  j  n(n  *  nPQ,M  .  Vl-Si?!1 

8Pl  1=1  1  2(»ls  -  1)  («LS  -  1) 

The  data  required  to  implement  this  model  are: 


Data  Requirement 

The  times  between  error  occurrences,  i.e.,  the  x^'s. 

4.3.3  Thompson  and  Chelson's  Bayesian  Reliability  Model 

The  last  model  considered  in  the  Bayesian  framework  for  software  reliability 
is  one  proposed  by  W.  E.  Thompson  and  P.  0.  Chelson.70  The  model  they  developed 
is  one  step  in  the  direction  of  obtaining  total  system  reliability.  Their  ulti¬ 
mate  goal  of  system  reliability  included  system  malfunctions  not  only  due  to 
software  but  to  hardware  and  unknown  or  ambiguious  source-related  malfunctions  as 
well.  In  a  paper  by  R.  Haynes  and  W.  E.  Thompson,71  this  total  system  reliability 
model  is  formulated.  The  one  aspect  of  this  model  that  this  paper  presents  is 
the  reliability  model  developed  for  the  software  related  errors.  This  model  at¬ 
tempts  to  account  for  the  fact  that  a  given  software  program  might  be  error-free 
(hence,  an  infinite  MTBF)72  and  it  provides  for  software  redesign  and  repair 
after  malfunctions  are  observed  in  a  given  test  phase.  The  specific  assumptions 
for  their  model  are: 
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Model  Assumptions 


(a)  The  program  is  not  corrected  during  a  testing  cycle-only  at  the  comple¬ 
tion  of  a  cycle  and  before  the  start  of  a  new  one. 

(b)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(c)  The  software  errors  are  assumed  to  occur  at  some  unknown  constant  rate, 
A.  The  total  number  of  errors  observed  in  a  testing  cycle  of  length  T  follows  a 
Poisson  distribution  with  parameter  AT;  i.e., 

-AT^-sf . 

f(f.|A)  =  6  f-^iJ  1  f.=0,l,...  (4.628) 

I  i 

(d)  If  p  denotes  the  probability  that  the  software  contains  one  or  more 
errors,  it  can  be  assumed  that  p  has  a  prior  distribution  that  is  beta  with 
parameters  a  +  1  and  b  +  1,  i.e., 

g(p)  =  — F-(-a  +  k  +  ^ —  pa(l  -  p)*3  0  <  p  <  1  a,  b  >  -1.  (4.629) 

T(a  +  l)r(b  +  1) 

The  parameter  a  is  thought  of  as  the  number  of  previously  delivered  software  pack¬ 
ages  with  errors  among  a  total  of  a  +  b  delivered. 

(e)  The  uncertainty  about  the  parameter  A  is  expressed  as  a  prior  distribu¬ 
tion  for  A.  It  is  assumed  gamma  with  parameters  To  and  fo+1;  i.e., 

T0(AT0)f° 

h(A)  =  -r^-  +  jy  exp  (-AT0)  A  >  0  .  (4.630) 

The  fo  can  be  thought  of  as  the  number  of  software-related  system  malfunctions  in 
previous  testing  of  total  duration  To. 

Thompson  and  Chelson  consider  two  situations.  One  is  the  situation  when 
it  is  known  before  testing  begins  that  the  software  contains  errors,  i.e.,  p  =  1. 
The  other  situation  is  when  there  is  uncertainty  about  whether  the  software  does 
or  does  not  contain  errors.  This  is  expressed  by  the  use  of  the  prior  in  assump¬ 
tion  (d).  If  ill  this  latter  situation,  an  error  is  discovered  in  the  testing 
cycle,  p  is  set  equal  to  1  and  the  prior  g(p)  is  made  a  Dirac  delta  function  at 

p  =  0. 

For  the  first  situation,  Thompson  and  Chelson  show  that  if  f^  errors  are 
observed  in  testing  time  T^,  then  the  posterior  distributions  for  A  and  R  (the 
software  reliability),  are: 


J 
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f.+f0+l  f.+fo 

(T  +  T0)  \ 

=  r(f.  +  f0  +  i) 


exp  [-A(T^  +  To)]  A  >  0  (A. 631) 


and 


f  (R|  f  s.)  = 


\  f o+£ i 

T.  +  T0  \ 

V  t  / 


An  (- 


0 


f0+f.-l  R 


/T.  +  T  \ 

’  V 

T  (f0  '+  TJ~~ 


(4.632) 


0  <  R  <  1. 


(The  t  is  th:  postulated  mission  time  for  the  program.)  The  distribution  for  R 
reflects  the  posterior  view  of  the  program  reliability  after  f^  errors  are  ob¬ 
served  in  the  current  testing  cycle. 


The  second  situation  is  the  one  in  which  no  errors  are  observed  in  the  test¬ 
ing  cycle  i;  i.e.,  f.  =  0.  This  generates  the  uncertainty  about  whether  the 
program  does  or  does  not  have  any  errors  still  residing  in  the  code.  For  this 
situation,  the  posterior  cumulative  distribution  functions  for  A  and  R  are  shown 
to  be: 


H(A|  f.  *  0)  = 


_ a  +  1 

a  +  b  +  2 


A  fo+1  fo  /  \ 

(T.  +  T0)  x  expf-x(Ti  +  T0)jdx 


(4.633) 


and 


F(R|f.  =  0,p)  =  p 


/ 


(T,  +  T0)f°  +  1 


r(f0  +  i) 


I  An  -  j 

L_  xj. 


T.+To-l 


d  x  0<R<1 


=  1 


R=1 . 


(4.634) 

(4.635) 


If  a  squared  error  loss  function  is  assumed  in  estimating  A  and  R,  the  Bayes 
estimates  for  A  and  R  are  then  the  means  of  the  respective  posterior  distribu¬ 
tions.  They  are  shown71  to  be: 


Estimates  -  Bayes 

l  -  (a  +  1)  .  (fo  +  1)„  (4.636) 

(a  +  b  +  2)  (T±  +  T0) 
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Data  Requirements 


(a)  The  number  of  software  errors  discovered  in  each  period  of  testing, 
i.e. ,  f . 's. 

(b)  The  length  of  testing  time  for  each  period,  i.e.,  the  T/s. 

(c)  For  the  total  number  of  software  packages  that  have  been  released,  the 
number  found  to  contain  errors.  These  numbers  are  used  in  determining  the  prior 
distribution  for  p  at  any  stage. 
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4.4  MARKOV  MODELS 

This  next  class  of  models  views  the  software  correction  process  as  a  dis¬ 
crete  space  system  in  which  a  transition  from  one  state  to  another  occurs  when¬ 
ever  an  error  detection  or  correction  is  made.  These  models  attempt  to  achieve  a 
more  precise  and  realistic  error  behavior  prediction  but  at  the  cost  of  a  great 
deal  of  added  complexity.  In  fact,  the  models  in  many  cases  cannot  be  used  to 
derive  a  closed  form  solution.  Only  large  sample  approximations  can  be  given  or 
approximate  numerical  solutions  can  be  stated.  Much  research  is  still  needed  in 
this  area  of  modeling.  This  section  is  included  for  completeness  in  the  presenta¬ 
tion  of  the  various  approaches  to  software  reliability  modeling.  Because  of  the 
complexity  of  the  models,  this  section  is  not  developed  in  great  detail.  The 
reader  is  referred  to  the  respective  research  articles  for  additional  details. 


4.4.1  Trivedi  and  Shopman's  Many  State  Markov  Models 

The  basic  model  and  its  generalizations  are  presented  in  a  paper  by  Ashok 
Trivedi  and  Martin  Shooman73  under  contract  to  the  Office  of  Naval  Research  and 
the  Rome  Air  Development  Center.  The  model  is  used  in  providing  estimates  of  the 
reliability  and  availability  of  a  software  program  based  upon  an  error  detection 
and  correction  process.  Availability  is  defined  as  the  probability  that  the  pro¬ 
gram  is  operational  at  a  specified  time.  The  software  can  be  viewed  in  either 
one  of  two  states  "up"  or  "down."  The  system  is  in  an  up  state  if  no  errors  have 
occurred  or  an  error  has  just  been  corrected.  The  software  is  in  a  down  state 
when  an  error  has  been  discovered  and  is  being  corrected.  The  sequence  of  up 
state  is  denoted  by  (n,  n  -  1,  n  -  2,...,n  -  k,...)  while  the  sequence  of  down 
states  is  denoted  as  (m,  m  -  l,...,m  -  k,...).  The  system  is  in  the  up  state, 
n  -  k,  if  the  (k  -  l)st  error  has  been  detected  but  the  kth  has  not.  It  is  in 
the  down  state,  m  -  k,  if  the  kth  error  has  been  detected  but  not  yet  corrected. 


The  specific  assumptions  for  the  model  are: 


Model  Assumptions 


(a)  The  transition  probability  from  state  i  to  state  j  (p^j)  is  dependent 

only  on  those  states  and  is  independent  of  all  past  states  except  the  last  one. 
(The  Markov  property.) 

(b)  The  error  detection  rate  for  the  state  n  -  k  is  known;  denote  it  as 

\  .  .  The  error  correction  rate  for  state  m  -  k  is  known;  denote  it  as  u  ,  . 
n-k  m-K 

(c)  Finite  nonzero  times  are  spent  by  the  system  only  in  the  system  states. 
The  transition  times  are  infinitesimally  small  so  the  probability  of  two  or  more 
error  detections  or  corrections  within  this  time  frame  is  zero. 


(d)  The  software  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 
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(e)  When  a  software  error  is  corrected,  it  is  done  without  the  introduction 
of  additional  errors. 

(f)  The  program  is  assumed  to  be  fairly  large  (the  order  of  10s  words  or 
more  of  code). 

Some  of  the  generalizations  of  this  basic  model  include  allowing  error  intro¬ 
duction  in  the  correction  process  and  the  system  can  be  in  more  than  two  states. 
A  generalization  of  the  basic  model  allows  a  third  state,  a  "noncritical  down" 
state.  The  reader  is  referred  to  Trivedi  and  Shooman's  paper  for  details. 

The  two  specific  cases  of  the  basic  model  that  this  pap.  reviews  are  called 

Model  I  and  Model  II.  For  Model  I,  the  error  detection  rate,  K  .  ,  and  the  error 

’  n-K 

correction  rate,  u  .  ,  are  taken  as  functions  of  the  number  of  errors  that  have 
’  rm-k’ 

occurred,  i.e.,  k.  In  Model  II,  the  rates  are  taken  as  functions  of  time.  The 
choice  between  the  two  models  is  determined  by  the  way  the  error  data  are  col¬ 
lected.  For  Model  I,  the  individual  errors  are  recorded  along  with  the  time  of 
occurrence  of  each  error.  For  Model  II,  the  number  of  errors  is  recorded  over  the 
operating  time  of  the  program. 

For  either  model,  the  derivation  of  the  reliability  and  availability  is  as 
follows.  Suppose  Pn_^(t)  denotes  the  probability  that  at  time  t  the  state  is 

in  the  n  -  k  up  state.  Similarly,  P  .  (t)  is  the  corresponding  probability,  the 

m-K 

system  is  in  the  m  -  k  down  state.  Then  the  availability  of  the  program  is: 

A(t)  =  P  {system  is  up  at  time  t}  (4.642) 

=  PQ(t)  +  P^jCt)  +  ...  (4.643) 

=  £  Pn  M  .  (4.644) 

k=0  n  K 

Thus,  only  the  probabilities  for  the  various  up  states  for  the  system  are  needed 
to  derive  the  availability.  The  reliability  on  the  other  hand,  depends  upon  the 
stage  of  debugging  since  the  smaller  the  number  of  residual  errors,  the  less 

likely  it  is  for  the  program  to  "discover"  them.  Suppose  the  system  has  just 
entered  the  state  n  -  k  at  time  t.  Suppose  this  time  is  renamed  as  t  =  0,  then 
in  the  interval  (0,  T^) ,  where  T^  is  the  time  of  discovery  of  the  kth  error,  the 

error  occurrence  rate,  \(k),  is  a  constant.  The  reliability  function  is  then 

R(t)  =  e'X(k)x  0  <  t  <  Tk,  k  =  1,2,....  (4.645) 

Hence,  after  the  (k  -  l)st  error  has  been  corrected,  only  \(k)  is  needed  to  esta¬ 
blish  the  reliability  for  the  program  for  all  times  between  the  occurrence  of  the 
(k  -  l)st  and  the  kth  error. 
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Now  consider  how  the  state  probabilities  are  derived.  First,  Model  I  is 
considered  with  the  special  case  of  constant  error  detection  and  correction  rates, 
i.e. , 


and 


Vk  =  Vk(k)  =  x 


*"W-k  =  Vk(k)  =  » 


k  =  0,1, 


k  =  0,1,. .. 


(4.646) 


(4.647) 


For  any  At  (At  small) ,  the  following  system  of  equations  represents  the  transi¬ 
tion  behavior  of  the  Markov  system: 


Pn(t  +  At)  =  (1  -  XAt)Pn(t), 


Pn_k(t  +  At)  =  (1  -  AAt)Pn.w(t)  +  pAt  Pm„t+1(t)  k  =  1,2,.. 


•n-k 


m-k+1 


(4.648) 

(4.649) 


and 


Pm-k(t  +  At)  =  (1  "  MAt)Pm-k(t)  +  XAtPn-k(t)  k  =  O*1’***  *  (^-650) 


By  dividing  both  sides  of  the  previous  equations  by  At  and  letting  At+0,  the  fol¬ 
lowing  set  of  differential  equations  is  obtained: 


Pn(t)  =  -\Pn(t). 


(4.651) 


Vk(t)  +  =  M  Vf.W 


and 


m-k(t) 


n-k 


m-k+1 


k  =  1,2,... 


+  (J  P  .  (t)  =  \  P  ,  (t) 
r  m-kv  ' 


n-k' 


k  =  0,1,2...  .  (4.652) 


Using  the  initial  conditions: 
Pu(0)  =  1  - 

W>  = 0 


k  =  1,2,3, . . 


and 


Vk'0’  = 0 


k  =  0,1,2,... 

Trivedi  and  Shooman73  show  that  the  solutions  to  this  system  of  equations  are: 

k-j 


P  .  (t)  = 
n-k 


.  \k  -\t  k  ,  .  k-j  f 


k+1  -(p-X)t 

c,  .e 
kj 


+  (-1)J  d.  . 

kji 


k  =  0,1,2, . . 
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the  rates  considered  as  functions  of  time,  tN;re  is  an  analogous  set  of  differen¬ 
tial  equations  as  obtained  in  Model  I;  namely, 


Pn(t)  =  -A(t)  Pn(t)  , 


Vk(t)  +  X(t)  Pn-k(t)  =  M(t)  Pm-k-H(t)  k  = 


and 


Pm-k(t)  +  |j(t)  Pm-k(t)  =  X(t)  Pn-k(t);  k  =  0>1*2 . 

The  initial  conditions  are  the  same  as  those  as  given  in  Model  I. 


(4.663) 

(4.664) 

(4.665) 


The  solutions  to  this  system  of  equations  may  not  even  exist  in  closed  form. 
As  for  Model  I,  numerical  solutions  must  therefore  be  relied  upon.  Trivedi  and 
Shooman's  paper  is  referred  to  for  additional  details. 

The  data  required  to  implement  this  model  are: 


Data  Requirements 

(a)  The  error  detection  rate  between  the  times  of  error  occurrence  which  is 
either  expressed  as  a  function  of  the  number  of  errors  detected  or  as  a  function 
of  time. 

(b)  The  error  correction  rate  between  the  times  of  error  occurrence  which  is 
expressed  either  as  a  function  of  the  numbers  of  errors  detected  or  time. 

The  interesting  aspect  of  this  model,  aside  from  the  Markovian  aspect,  is 
that  no  parameters  are  estimated.  The  error  detection  and  correction  rates  are 
needed  as  input  into  the  model  formulation.  If  these  are  unknown,  which  is  the 
usual  situation,  they  need  to  be  estimated.  In  the  example  application  considered 
by  Trivedi  and  Shooman,  empirical  estimates  of  the  rates  are  obtained  as  p  func¬ 
tion  of  time  by  using  the  number  of  software  error  reports  per  month  for  the  de¬ 
tection  rate  and  the  number  of  closed  software  error  reports  per  month  for  the 
correction  rate.  Sukert17  uses  the  number  of  errors  found  and  corrected  per  day 
as  reported  in  software  error  forms  to  estimate  the  rates  in  this  application.  It 
employs  standard  regression  analysis  as  well  as  fitting  some  nonlinear  functions 
to  the  data  to  estimate  the  curves  \(t)  and  m(t). 


4.4.2  Littlewood's  Semi-Markov  Model 


The  last  model  considered  in  this  paragraph  was  proposed  by  Bev  Little- 
wood.74’  75  The  model  incorporates  the  structure  of  the  program  in  developing  its 
availability.  One  of  the  major  weaknesses  of  the  previous  time-dependent  models 
is  that  the  structure  of  the  program  is  not  considered  in  determining  its  reli¬ 
ability.  Littlewood  adopts  a  modular  approach  to  the  software  and  att^ipts  to 
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describe  this  structure  via  the  program's  dynamic  behavior  using  a  Markov  assump¬ 
tion.  The  program  comprises  a  finite  number  of  modules  with  exchanges  of  control 
between  them  that  follow  a  semi-Markov  law.  The  time  spent  in  a  given  module  can 
be  taken  as  a  random  variable  with  any  distribution  (hence,  semi-Markov)  which  is 
characteristic  of  the  module  and  the  module  that  it  transitions  to.  The  specific 
assumptions  of  this  model  are: 


Model  Assumptions 

(a)  The  program  is  composed  of  M  modules.  Transitions  between  modules  are 
such  that  the  probability  that  the  program  terminates  one  module  to  enter  another 
is  independent  of  the  time  the  first  module  is  entered  (semi-Markov  property). 


(b)  When  the  program  is  in  module  i,  the  failures  are  assumed  to  follow  a 
Poisson  process  with  parameter  v^. 

(c)  When  module  i  calls  module  j  with  probability  p..,  the  probability  of  a 
failure  in  the  interface  between  the  two  modules  is  cl  j . 

(d)  The  distribution  of  the  time  sp^nt  in  module  i  before  entering  module  j 
depends  upon  only  i  and  j  and  is  known  only  via  the  first  two  moments  (Jj1"*  and 
M2ij. 


(e)  The  program  is  operated  in  a  similar  manner  as  the  anticipated  opera¬ 
tional  usage. 

(f)  Each  failure  results  in  a  random  variable  cost.  The  random  variables 
are  assumed  independent  with  distributions  dependent  on  the  module  or  interface  in 
which  the  failure  occurs.  The  distributions  are  only  known  through  their  first 
two  moments. 

The  last  assumption  is  optional;  it  is  only  needed  if  an  overall  failure  cost 
analysis  is  needed. 

Suppose  N(t)  is  the  total  number  of  failures  (both  within  modules  and  be¬ 
tween)  observed  in  the  program  in  the  time  interval  (0,t).  Deriving  the  distri¬ 
bution  of  N(t)  for  a  specified  t  is  an  extremely  difficult  if  not  impossible  task. 
A  complete  description  of  the  behavior  of  N(t)  requires  knowledge  of  the  distribu¬ 
tions  of  times  within  modules,  a  requirement  that  is  usually  unattainable  in  prac¬ 
tice.  Littlewood  derives  an  asymptotic  result  pertaining  to  the  behavior  of  N(t). 
If  the  very  plausible  assumption  for  a  modular  program  is  made  that  the  individual 
failure  rates  are  much  smaller  than  the  switching  rates  between  modules,  then  the 
failure  point  process  of  the  integrated  program  is  asymptotically  a  Poisson  pro¬ 
cess  with  rate  parameter 
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M  M 


.E  S  n.p. .  (pi^v.  +ot..) 

,£LJ7.1...iZ?1.,:.H!_  1 - Hi  (4.666) 

MM.. 

£  £ 

1=1  1=1  j 

as  v.,  a..  ->0  where  n  =  (II.)  satisfies  II  •  P  =  n  with  52  n.  =  1  and  P  is  the 

i  ’  ij  ~  4  i*  ~  ~  fpf.  i 

MxM  transition  matrix  of  the  system. 

The  interesting  aspect  of  this  result  is  demonstrated  by  rewriting  the  pre¬ 
vious  expression;  i.e.,  if 


a  «  _ ifl  _  J _ 

i  MM  ,  . 

E  E  n.p^Mi1*1 

i=l  j=l  1 


(4.667) 


«B*  (5  "-I 
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b.  .  =, - 

1J  M  M 


nipij 


(4.668) 


E  E  nipijMi1J 

i=l  j=l  1 

then  the  previous  expression  can  be  reexpressed  as: 
M  MM 

E  a4v,*  +  £  £  b..a. .  . 

i=l  1  1  i=l  i=l 


(4.669) 


The  a^  represents  the  limiting  proportion  of  time  spent  in  module  i  while  b^  is 

the  limiting  frequency  of  i  to  j  module  transfers.  It  is  often  possible  to  esti¬ 
mate  them  directly.  An  extremely  complex  description  of  the  behavior  of  N(t)  is 
therefore  represented  asymptotically  in  a  very  simplistic  manner. 

Suppose  interest  is  also  in  a  failure  cost  analysis.  If  Y^(t)  represents  the 

random  variable  for  the  cost  of  a  failure  in  module  i,  and  (t)  represents  the 

random  variable  cost  of  a  failure  in  transfer  from  module  i  to  j,  the  total  pro¬ 
gram  cost  is: 


M  M  M 

Y(t)  =  £  Y.(t)  +  £  £  Y..(t)  . 

i=l  1  i=l  j=l 


(4.670) 
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Again  the  exact  description  of  Y(t)  is  extremely  complex  if  not  impossible  to  de¬ 
velop.  Littlewood  develops  an  asymptotic  result  which  depends  only  upon  the  means 
and  variances  of  the  defining  distributions.  He  shows  that: 

Y(t)  -  pt 


Of 


* 


N(0,1) 


(4.671) 


where  |j  is  the  (asymptotic)  mean  cost  incurred  per  unit  time  for  the  integrated 
program  and  is: 


M  M 

_  £  Miij) 

MM 

E  E  nipiiMr 

i=l  j=l  1  1J 


(4.672) 


ij 


>5 

h 


where  p*1  and  pi1*-*  are  the  means  of  Y^,  Y^  respectively. 

The  reader  is  referred  to  Littlewood 's  papers  for  additional  details  pertain¬ 
ing  to  this  result  and  a  definition  of  <J. 

The  major  problem  with  this  model,  as  with  the  previous  case  is  that  all  of 
the  parameters  that  make  up  the  model  are  input;  therefore,  they  must  be  known  or 
estimated.  The  data  required  for  this  model  are: 


Data  Requirements 

(a)  The  transition  probabilities  from  modules  i  to  j ,  i.e.,  the  p^'s. 

(b)  The  error  rates  within  the  modules  i,  i.e.,  the  v^'s. 

(c)  The  first  two  moments  of  the  distribution  for  the  time  spent  in  module 
i  before  transferring  to  module  j,  i.e,  the  Pi^'s  and  M21J*8* 

(d)  The  probabilities  of  failures  occurring  at  the  interfaces  between 
modules ,  i.e.,  the  a . . * s . 

(e)  If  a  cost  analysis  is  desired,  the  first  two  moments  of  the  cost  distri¬ 
butions  of  failure  within  and  between  modules. 


As  can  be  seen,  the  data  required  can  easily  prohibit  the  use  of  this  model. 
An  additional  factor  to  cousider  before  applying  this  model  (and  the  previous  one) 
is  assuring  that  the  Markov  property  is  satisfied.  This  could  prove  the  most 
formidable  problem  of  all  in  applying  Markov  Models. 
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CHAPTER  5 

COMPARISON  OF  RELIABILITY  MODELS 


For  this  part  of  the  report,  some  studies  and  the  results  pertaining  to 
comparing  the  performance  of  the  various  reliability  approaches  and  models  on 
software  error  data  sets  are  described.  Any  common  conclusions  that  have  been 
reached  among  the  studies  are  pointed  out.  The  studies  described  all  involve  the 
comparison  of  at  least  three  or  more  software  reliability  models.  Excluded  are 
various  individual  studies  that  have  been  done  on  a  given  software  model  only. 
This  includes  the  results  of  applying  Musa's  Model  to  error  data  sets  54,55  aud 
Littlewood's  Bayesian  Mode].  Musa's  Model  has  been  applied  to  actual  software 
error  data,  with  some  success.  Littlewood  also  has  applied  his  Bayesian  Model  to 
some  data  sets  with  a  good  match  obtained  between  the  predicted  and  actual  error 
observations.  However,  this  comparison  is  strictly  concerned  with  results  that 
can  be  given  in  the  performance  of  one  reliability  model  in  contrast  to  another. 
From  these  studies,  some  guidelines  can  be  established  for  employing  a  model  in-  a 
given  situation.  In  this  paragraph,  the  analysis  is  limited  to  the  Time  Domain 
Approach  rather  than  Error  Seeding/Tagging  and  the  Data  Domain  Approaches  since 
no  major  studies  have  been  done  to  compare  the  performance  of  the  different  ap¬ 
proaches.  The  first  major  effort  was  made  by  Alan  Sukert  of  the  Rome  Air  Develop¬ 
ment  Center. 17,76,77  The  study  involves  five  major  models:  The  Jelinski-Moranda 
Model  (Paragraph  4.2.3),  the  Schick-Wolverton  Model  (Paragraph  4.2.4),  a  modified 
Schick-Wolverton  Model  (Paragraph  4.2.4. 1),  and  Geometric  (Paragraph  4.2.6)  and 
Modified  Geometric  Models  (Paragraph  4.2.6. 1).  Other  models  were  considered  but, 
due  to  data  requirements  for  these  models  (e.g.,  CPU  time),  could  not  be  used. 
The  chosen  models  were  applied  to  four  large  scale  DOD  software  projects.  Using 
the  software  error  reports  filed  for  each  of  the  projects,  the  error  counts  per 
day  and  per  week  were  used  as  input  into  the  software  models.  (All  models  were 
modified  to  allow  more  than  one  error  per  time  frame.)  The  study  considered 
estimates  for  the  respective  models  obtained  from  maximum  likelihood  and  least 
square  procedures. 

The  basic  conclusions  drawn  from  this  comparative  study  are: 

(a)  The  grouping  by  weeks  does  better  in  predictive  ability  than  by  the  day. 

(b)  The  Jelinski-Moranda  and  Schick-Wolverton  Models  give  reasonable  predic¬ 
tions  for  small  projects  while  the  modified  Schick-Wolverton  Model  does  better  for 
larger  ones. 

(c)  The  Geometric  Models  are  better  to  use  when  the  MTBF  or  reliability 
estimates  are  of  concern. 

A  major  problem  experienced  in  the  application  of  these  models  in  this  study  (as 
in  the  others)  is  the  problem  of  convergence.  The  estimation  procedures  failed 
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in  many  instances  to  come  up  with  model  parameter  estimates.  This  is  discussed 
in  some  detail  after  the  next  study  is  described. 

The  second  major  comparative  effort  was  undertaken  by  Hughes  Aircraft  Com¬ 
pany.33  The  study  involves  the  Generalized  Poisson  Model  of  Paragraph  4.2.5  with 
g(Xi,...,Sn)  =  X® ;  (a  unknown),  a  binomial  model  and  the  Nonhomogeneous  Poisson 

Process  Model  of  Paragraph  4.2.9.  The  various  models  employ  both  maximum  likeli¬ 
hood  and  least  squares  in  the  estimation  of  model  parameters.  These  models  are 
applied  to  16  sets  of  electronic's  system  computer  program  software  data.  The 
major  conclusions  reached  are: 

(a)  Generally,  the  model  fits  to  the  data  are  poor,  but  the  best  fitting  of 
the  three  and  the  one  applicable  to  the  most  data  sets  is  the  Generalized  Poisson 
Model. 

(b)  Grouping  the  error  data  by  a  time  period  has  better  convergence  proper¬ 
ties  than  ungrouped. 

(c)  Maximum  likelihood  and  least  squares  estimates  for  a  given  model  are 
similar. 

The  major  problem  of  lack  of  convergence  to  parameter  estimates  was  experi¬ 
enced  in  this  study  as  well.  The  authors  suspect,  as  does  Sukert  in  the  previous 
study,  that  a  major  problem  causing  this  lack  of  convergence  is  violation  of  the 
assumptions  on  which  the  models  are  based,  especially  the  violation  of  a  nonin¬ 
creasing  error  rate.  The  Hughes  report  finds  that  by  plotting  the  estimated  error 
rate  whenever  it  is  increasing  in  a  region  is  precisely  the  region  in  which  con¬ 
vergence  problems  are  experienced.  The  types  of  convergence  problems  encountered 
include  lack  of  convergence,  oscillation,  convergence  to  a  nonoptimal  solution, 
and  nonuniqueness  of  the  solution.  These  problems  are  especially  experienced  by 
the  MLEs.  The  report  employs  a  second  derivative  criterion  to  weed  out  nonoptimal 
solutions,  but  the  report  points  out  that  this  cannot  b$  relied  on  completely 
because  of  computer  precision  problems  in  finding  th“  optimal  solutions. 

The  third  study  was  undertaken  by  Dayton  University  under  contract  to  the  Air 
Development  Center.18  The  study  applied  the  Jelinski-Moranda  Model,  the  Geometric 
Model,  the  extended  Jelinski-Moranda  Model  (Paragraph  4. 2. 3. 2),  and  Schneidewind* s 
Model  with  approach  (c)  (Paragraph  4.2.8)  to  software  error  data.  The  first  two 
are  applied  to  two  data  sets  in  which  the  times  between  error  occurrences  are 
recorded  while  the  latter  two  are  applied  to  two  data  sets  in  which  the  error 
counts  are  recorded.  The  conclusions  are: 

(a)  If  the  times  between  error  occurrences  are  available,  the  Geometric 
Model  does  a  better  predictive  job,  but  if  the  error  counts  per  time  interval  are 
available,  the  Schneidewind  Model  is  preferred. 

(b)  The  extended  Jelinski-Moranda  Model  and  Schneidewind' s  Model  give  simi¬ 
lar  results,  but  the  extended  Jelinski-Moranda  Model  is  very  sensitive  to  changes 
in  the  data. 
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(c)  The  Jelinski-Moranda  Model  tends  to  estimate  a  smaller  number  of  errors 
remaining  than  the  Geometric  Model,  illustrating  the  "optimistic"  tendency  of  the 
exponential  class  of  distributions. 

The  fourth  study  was  performed  by  the  University  of  Utah.32  In  contrast  to 
the  previous  three  studies,  this  study  compares  a  number  of  models  using  deter¬ 
ministically  generated  error  data  rather  than  actual  data.  The  times  of  error 
occurrences  are  generated  on  a  computer,  following  the  underlying  model  with  known 
model  parameters.  The  models  considered  include:  the  Jelinski-Moranda  Model,  the 
Geometric  Model,  and  Musa's  Model.  A  Monte-Carlo  study  of  the  behavior  of  the 
least  squares  and  MLEs  was  undertaken.  The  results  of  this  study  are: 

(a)  A  strong  positive  correlation  is  indicated  among  the  various  estimates 
for  total  number  of  errors, 

(b)  The  estimate  of  the  MTBF  is  best  for  the  Geometric  Model, 

(c)  The  accuracy  of  the  estimates  increases  as  either  the  total  number  of 
errors  increase  or  the  number  of  errors  remaining  decreases,  and 

(d)  The  least  squares  estimation,  using  the  times  between  error  occurrences 
rather  than  the  actual  times  of  the  occurrences,  does  not  perform  as  well  as  the 
other  estimators. 

The  last  few  paragraphs  summarize  the  results  of  the  various  studies  under¬ 
taken  to  compare  the  various  models.  As  can  be  seen  from  the  studies,  additional 
comparative  research  is  needed.  Many  of  these,  studies  employed  error  data  that 
were  gathered  without  the  data  requirements  or  assumptions  of  the  various  models 
in  mind.  What  is  needed  is  a  large  scale  effort  in  which  the  data  are  gathered 
under  a  controlled  environment.  Currently  such  a  study  involving  the  Nonhomo- 
genous  Poisson,  the  IBM  Poisson  Model,  the  Generalized  Poisson  Model,  the 
Jelinski-Mcranada  Model,  and  the  Geometric  Poisson  Model  is  being  undertaken 
by  Hughes  Aircraft  for  the  Rome  Air  Development  Center.  The  data  are  being  col¬ 
lected  specifically  for  software  reliability  applications.  Although  the  final 
report  is  not  written  yet,  an  interim  report  (Reference  73)  finds  some  results 
similar  to  the  previous  studies.  Specifically,  the  major  problem  of  convergence 
and  the  violation  of  the  model  assumptions  are  found.  Again  the  major  violation 
is  an  error  rate  that  is  nonconstant  during  a  testing  interval  and  nondecreasing 
over  all  intervals. 

It  is  difficult,  based  upon  the  results  of  the  studies,  to  provide  clear 
cut  guidelines  in  applying  the  software  models.  We  can  only  conjecture.  Out 
of  the  various  models  considered,  it  appears  that  the  Generalized  Poisson  or 
Schneidewind' s  Model  approach  (c)  might  be  best  suited  for  count  data.  The 
Geometric  Models  should  be  considered  when  estimating  MTBF.  The  convergence 
problem  appears  to  diminish  as  the  length  of  the  testing  period  increases.  How¬ 
ever  there  is  no  method  to  determine  what  the  optimal  length  of  a  testing  in¬ 
terval  should  be.  This  depends  upon  the  underlying  error  generation  process  which 
is  not  known. 
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Many  of  the  models  discussed  in  this  report  have  yet  to  be  compared  on  the 
basis  of  performance  with  others.  None  of  the  Bayesian  Models,  Markov  Models, 
the  Error  Seeding/Tagging  Models,  or  the  Data  Domain  Approach  to  software  reli¬ 
ability  modeling  have  been  included  in  any  comprehensive  study.  (Note:  Hughes' 
current  research  will  incorporate  Lit.tlewood' s  Bayesian  Model).  This  report  can 
only  present  what  has  been  done  and  what  those  limited  results  indicate.  Much  is 
yet  to  be  done-if  it  even  can  be  done.  A  large  scale  controlled-data  collection, 
in  which  the  CPU  time  and  wall  clock  time  are  simultaneously  gathered  for  the 
purpose  of  comparing  as  many  different  models  as  possible,  may  be  economically 
and  administratively  infeasible.  Moreover,  for  the  modeling  of  software  error 
generation,  no  one  model  is  applicable  in  all  instances.  The  software  analyst 
needs  a  collection  of  software  models  which  have  demonstrated  themselves  in 
various  environments  and  comparative  studies.  From  this  collection,  the  analyst 
judiciously  selects  the  one  most  applicable  to  his/her  situation.  Flexibility  and 
adaptability  are  the  keys  to  successful  modeling. 
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CHAPTER  6 


"QUICK”  ESTIMATES  OF  SOFTWARE  RELIABILITY  MEASURES 


This  last  paragraph  briefly  presents  some  proposed  "quick"  estimates  for 
various  software  reliability  measures.  These  procedures  do  not  require  the  exten¬ 
sive  error  data  base  of  the  previous  sections.  The  view  taken  toward  software 
reliability  is  very  simplistic  and  pragmatic  in  nature.  The  two  procedures  dis¬ 
cussed  are  not  advocated  in  this  report,  but  are  included  for  completeness.  The 
purpose  of  this  report  is  to  review  all  of  the  various  procedures  that  are  advo¬ 
cated  in  determining  the  reliability  status  of  a  set  of  software. 

6.1  MTBF  ESTIMATION 

This  very  simple  measure  of  MTBF  was  proposed  by  Gregory  Hansen  of  Systems 
Engineering  Laboratories.79  When  a  software  program  is  first  released,  there  are 
only  a  few  users  and  hence  the  failure  generation  is  a  minimum.  This  means  that 
the  MTBF  is  fictitiously  high,  giving  the  software  manager  a  false  sense  of  secur¬ 
ity.  As.  the  software  begins  to  be  used,  the  MTBF  can  be  expected  to  rise  slightly 
as  the  initial  gross  errors  are  discovered  and  eliminated.  However,  in  later 
years  as  more  and  more  users  test  the  software,  the  MTBF  drops  significantly. 
Finally,  the  software  reaches  a  "mature"  state  and  the  MTBF  increases  sharply. 

This  behavior  is  reflected  in  the  following  formula: 


(H.  +  N.  ,  +  ...  +  N.  ..,)*C 
1  1-1  1-1+1 


(6.1) 
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where  M^  =  number  of  software  errors  discovered  in  year  i, 

MTBF(i)  a  MTBF  for  year  i, 

=  Number  of  copies  of  the  program  in  use  for  year  i, 

I  =  Number  of  years  the  program  has  been  used, 

and  C  =  Estimates  of  the  average  number  of  hours  that  the  product  is  used 

in  a  year. 

An  example  calculation  is  as  follows.  Suppose  the  MTBF  is  desired  for  year  i  =  3 
with  the  software  being  distributed  to  users  for  2  years.  During  the  third  year, 
there  has  been  a  total  of  10  users  and  during  the  second,  a  total  of  5  users.  The 
estimated  average  number  of  hours  that  the  program  is  used  during  the  third  year 
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is  500  hours.  During  the  third  year,  a  total  of  100  errors  have  been  observed. 
The  formula 


(N3  +  N2)*C 
“O)  ",  - — 


(10  +  5)*500 
100 


75  hrs 


(6.2) 


reflects  the  total  estimated  number  of  hours  of  use  by  all  users  during  the  given 
year,  divided  by  the  number  of  errors  found  during  the  year.  The  function  can  be 
plotted  against  time  to  see  the  error  behavior  of  the  software  package.  Once  the 
program  reaches  maturity,  future  values  of  MTBF  can  be  predicted  by  extending  the 
curve. 


6.2  PRAGMATIC  SOFTWARE  RELIABILITY  ESTIMATION 

The  last  method  considered  is  proposed  by  John  Wall  and  Paul  Ferguson.80 
Using  the  basic  premise  that  the  failure  rate  of  software  decreases  as  more  soft¬ 
ware  is  used  and  tested,  they  formulate  a  relationship  between  the  number  of  fail¬ 
ures  and  the  ’'maturity"  of  the  software.  Specifically,  the  relationship  proposed 
is: 

c  =  Co  (y“  (6.3) 

where  C  is  the  cumulative  number  of  errors  experienced  for  a  software  program  of 
maturity  M.  Co  and  a  are  constants  determined  empirically  by  plotting  the  cumu¬ 
lative  number  of  errors  versus  the  maturity  level  of  the  software.  Mo  is  a  scal¬ 
ing  constant.  Typically,  the  units  of  M  and  Mo  are  expressed  as:  amount  of 
calendar  time  expended,  processor  or  QPU  time,  man-months  of  testing,  or  the  num¬ 
ber  of  tests  executed. 


The  failure  rate,  R,  is  then  determined  as 


R 


(g- 


a-1 


(6.4) 


For  convenience  this  is  expressed  as: 


(6.5) 


where  Ro  is  simply  a  constant.  Again  the  terms  Ro  and  o<  can  be  determined  empiri¬ 
cally  from  the  data.  For  example,  failures  per  CPU  second  can  be  plotted  versus 
number  of  CPU  seconds  of  operation  to  determine  Ro  and  a.  Care  must  be  taken  to 
ensure  consistency  of  the  units  in  the  functional  relationships. 

The  application  of  this  method  is  applied  to  a  number  of  data  sets  in  their 
paper.  The  reader  is  referred  to  that  paper  for  additional  details. 
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CHAPTER  7 

SUMMARY  AND  CONCLUSIONS 


With  the  ever  increasing  role  that  software  is  playing  in  the  weapon  systems 
and  the  increased  complexity  of  the  programs  because  of  that  role,  a  dramatic  in¬ 
crease  in  the  cost  of  the  software  over  the  life  cycle  of  the  weapon  system  is 
seen.  Greater  emphasis  has  thus  been  placed  in  determining  more  cost  effective 
ways  of  software  development  and  testing.  One  such  method  that  has  developed  over 
the  last  10  years  is  the  calculation  of  a  software's  reliability.  By  having  a 
quantative  measure  of  a  program's  or  a  program  module's  reliability,  a  software 
manager  can  best  determine  the  allocation  of  testing  personnel  and  just  how  much 
testing  to  employ  before  release  to  the  user. 

This  report  provides  a  review  of  the  various  approaches  to  estimating  that 
reliability.  The  three  major  approaches  are  categorized  as:  Error  Seeding/Tagg¬ 
ing,  Data  Domain,  and  Time  Domain.  The  Error  Seeding/Tagging  Approach  u  es  the 
concept  of  error  introduction  into  the  software.  Based  upon  the  number  of  in¬ 
serted  errors  and  inherent  errors  found  in  the  testing  phase,  the  total  number 
of  errors  still  residing  in  the  program  can  be  estimated.  The  major  problem 
with  this  approach  is  the  implementation.  How  are  errors  of  the  same  nature  and 
distribution  as  the  inherent  errors  inserted  into  a  program?  The  Data  Domain 
Approach  bases  the  reliability  estimation  on  the  number  of  successful  execution 
runs  out  of  the  total  number  of  runs  attempted.  In  addition,  the  approach  tries 
to  incorporate  the  input  domain  structure  into  the  estimation  process.  The  input 
space  is  broken  down  into  regions  which  are  assigned  probabilities  based  upon 
anticipated  operational  profiles.  Random  samples  from  the  input  space  are,  then 
drawn  according  to  these  probabilities  and  the  count  of  successful  runs  made 
from  them  . are  used  in  the  reliability  calculation.  The  major  weakness  with  this 
approach  is  the  stratification  of  the  input  space  and  the  resulting  probability 
assignments . 

The  last  approach,  which  this  paper  deals  with  the  most,  is  the  Time  Domain. 
This  approach  attempts  to  model  the  error  generation  process  as  observed  over  time 
(either  CPU  or  wall  clock) .  This  is  done  using  the  time  of  error  occurrence  (or 
equivalently,  the  time  between)  or  the  number  of  errors  observed  over  a  testing 
interval.  Many  of  the  models  are  based  upon  an  underlying  Poisson  process  for  the 
error  generation  over  a  specified  time  frame  or  an  exponentially  distributed  ran¬ 
dom  variable  for  the  time  between  error  occurrences.  The  Time  Domain  Approach  can 
itself  be  categorized  into  three  types  of  models:  "Classical,"  "Bayesian,"  and 
"Markov."  The  "Classical"  Models  can  be  traced  back  to  their  origin  within  hard¬ 
ware  reliability  theory.  Many  of  the  concepts  of  hardware  reliability  theory 
(MTBF,  hazard  rate,  reliability  function)  are  adapted  to  the  field  of  software. 
Moreover,  models  of  this  class  tend  to  view  the  errors  inherent  in  a  program  to  be 
of  the  same  order  of  magnitude  and  the  correction  of  any  one  of  them  has  the  same 
order  of  impact  upon  the  program.  The  "Bayesian"  viewpoint  takes  this  impact  and 
treats  it  as  a  random  variable.  It  is  not  known  what  effect  the  correction  of  an 


error  might  have  upon  the  behavior  of  the  program.  When  errors  are  discovered 
early  in  the  testing  cycle,  it  is  expected  that  the  most  dramatic  improvement 
in  the  performance  of  the  program  occurs  after  their  correction.  Errors  discover¬ 
ed  late  in  the  cycle  have  the  least  dramatic  improvement.  The  "Markov"  Models 
attempt  to  formulate  the  error  generation  process  over  time  as  a  Markov  process 
in  which  transition  probabilities  are  either  given  or  derived.  These  probabili¬ 
ties  are  the  state  transition  probabilities  for  moving  from  one  state  to  another. 

As  in  the  previous  approaches,  the  Time  Domain  also  has  its  share  of  prob¬ 
lems.  The  Markov  Models  are  extremely  complex  and  difficult  to  apply.  Most 
of  the  results  are  either  for  special  cases  or  are  asymptotic  in  nature.  The 
Bayesian  Models  represent  a  more  realistic  approach  to  modeling  the  actual  error 
generation/ correction  process.  The  difficulty  here,  as  with  the  Bayesian  theory 
in  general,  is  the  specification  of  a  prior  distribution  for  the  error  rate.  In 
addition,  little  has  been  done  in  comparing  this  class  of  models  to  models  of  a 
"classical"  nature.  The  major  weakness  for  the  Classical  Models  is  an  oversen¬ 
sitivity  to  the  violations  of  the  assumptions  upon  which  they  rest.  They  are 
especially  sensitive  to  an  increasing  failure  rate  within  the  data.  This  in¬ 
creasing  rate  may  be  due  to  many  reasons;  introduction  of  new  errors  in  the 
correction  process,  nonuniform  testing,  and  nonuniform  application  of  testing  man¬ 
power  throughout  the  testing  cycle.  The  last  two  are  especially  common  occur¬ 
rences  in  typical  software  testing  programs .  Another  problem  that  these  models 
face  is  a  lack  of  independence  among  the  errors.  In  many  instances,  the  dis¬ 
covery  of  one  error  quickly  leads  to  others  generating  a  "clumping"  effect  of  the 
errors  over  time.  These  various  violations  lead  to  poor  fits  of  the  models  and 
convergence  problems  in  the  estimation  process. 

Various  studies  have  been  undertaken  to  compare  the  performance  behavior 
among  the  models,  but  no  clear  superior  model  has  arisen.  It  is  felt  by  this 
author  that  no  one  model  can  be  advocated  for  all  applications.  A  collection  of 
models  that  have  demonstrated  themselves  over  a  large  class  of  problems  should  be 
considered.  The  software  analyst  should  then  pick  from  this  class  the  one  that 
is  most  effective  in  modeling  his/her  set  of  data.  Modeling  has  always  been  an 
interative  procedure,  it  includes  choosing  a  candidate  model,  estimating  the 
parameters  of  the  model,  testing  the  adequacy  of  the  model,  and  cycling  back  if 
necessary. 

Much  research  is  yet  to  be  done  in  this  new  field;  however,  software  reli¬ 
ability  modeling  can  provide  an  effective  aid  to  the  software  manager.  Some  of 
the  applications  of  these  models  demonstrate  this.  What  is  to  be  kept  in  mind, 
however,  is  that  it  is  one  of  many  tools  available  in  developing  cost  effective 
software.  By  careful  consideration  of  the  collection  of  data  for  these  models,  to 
ensure  the  model  assumptions  are  satisfied  as  much  as  possible,  and  by  using  a 
collection  of  models  that  appear  "robust"  to  violations  of  assumptions  which  can¬ 
not  be  met,  the  models  provide  a  useful  aid.  Using  ths  common  techniques  of 
modeling,  the  chosen  model  can  be  a  useful  quantitative  measure  to  determine  the 
length  of  testing  and  manpower  utilization.  Otherwise,  if  the  manager  is  asked, 
"What  led  you  to  make  the  decision  to  release  the  software?"  What  can  he  say? 
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